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Abstract 

This  paper  develops  a  beam  pattern  design  procedure  for  general  multidimensional 
irregular  sonar  arrays  that  incorporates  the  not  well  understood  effects  of  array  geom- 
etry into  the  design  process.  The  procedure  is  implemented  by  generating  a  "penalty 
function"  in  a  spectral  covariance  function  form.  Processing  the  penalty  function  causes 
beam  pattern  high  sidelobes  to  be  penalized  and  the  main  lobe  to  be  emphasized.  This 
is  accomplished  by  forming  the  penalty  function  in  terms  of  an  isotropic  noise  field  of 
specified  strength  modified  with  a  finite  sector  of  low  coherent  energy  and  stabilized  with 
incoherent  sensor  noise.  By  inputting  the  penalty  function  into  a  minimum  variance 
beamformer,  the  beam  pattern  and  aperture  weights  are  calculated  based  on  the  given 
array  geometry.  The  beamformer  used  is  Capon's  Maximum  Likelihood  Method.  The 
array  used  to  test  the  procedure  is  located  on  a  sixty  degree  sector  of  a  cylindrical  sur- 
face. The  procedure  is  implemented  by  two  different  methods,  each  with  some  desirable 
characteristics.  One  method  suppresses  sidelobes  directly  by  the  placement  of  nulls.  The 
other  method  suppresses  sidelobes  indirectly  by  the  enhancement  of  the  main  lobe  with 
anti-nulls.  Both  methods  are  evaluated  in  terms  of  a  sensitivity  factor  which  constrains 
the  maximum  white  noise  array  gain.  Results  show  that  both  methods  result  in  sidelobe 
levels  that  range  from  20  to  35  dB  lower  compared  to  a  conventional  beam  pattern  with 
uniform  aperture  weighting  and  that  the  design  procedure  is  applicable  to  beam  patterns 
steered  to  both  true  broadside  and  to  off-broadside  directions. 
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Massachusetts  Institute  of  Technology 


Acknowledgements 

My  two  years  as  a  graduate  student  have  been  an  enjoyable  and  rewarding  period. 
Many  people  have  aided  me  in  the  course  of  my  studies,  more  than  I  can  reasonably 
mention  here.  However,  there  are  a  few  people  who  deserve  special  mention. 

First  and  foremost,  I  must  thank  my  advisor,  Dr.  Arthur  B.  Baggeroer.  Without 
his  invaluable  guidance,  encouragement  and  patience,  I  would  not  have  been  able  to 
successfully  complete  my  studies.  I  must  also  thank  Art  for  his  concern  and  interest  in 
my  career. 

A  special  thanks  also  goes  to  Eddie  Scheer,  computer  wizard  and  musical  guru.  Eddie 
had  to  contend  with  numerous  computer  related  questions  and  requests  while  maintaining, 
for  the  most  part,  his  jovial  attitude.  In  particular,  Eddie  was  an  available  outlet  when  a 
break  was  needed  and  was  an  immense  help  when  he  typed  a  portion  of  this  paper. 

I  am  also  indebted  to  Rob  Fricke.  Rob  was  the  one  person  who  was  always  available 
and  always  willing  to  discuss  my  studies  and  unrelated  topics. 

I  must  also  recognize  my  fellow  Joint  Program  adventurers.  Although  many  of  us  had 
(and  still  have)  the  same  goals,  we  have  achieved  those  goals  through  different  avenues 
and  in  doing  so,  have  acquired  a  shared  base  of  knowledge  not  found  anywhere  else. 

I  am  most  grateful  to  the  United  States  Navy  in  conjunction  with  the  Oceanographer 
of  the  Navy  for  allowing  me  this  valuable  opportunity  and  for  provided  the  financial 
backing. 

And  last,  I  must  recognize  the  love  and  support  I  have  received  from  my  family.  My 
wife,  Laurie,  and  my  daughters,  Brittany  and  Brianna,  have  kept  me  in  touch  with  reality 
and  my  true  priorities  in  life. 


Contents 


1  Introduction  9 

1.1  Background 9 

1.2  Objectives 10 

1.3  Overview 11 

2  Problem  Development  13 

2.1  Array  Geometry     13 

2.2  Array  Processing  Approach 15 

2.2.1  Signal  Model 15 

2.2.2  Spectral  Covariance  Matrix 17 

2.2.3  MLM  Development 21 

2.2.4  Array  Gain 24 

2.2.5  Spatial  Scanning 25 

2.3  Beam  Pattern  Design  Method 29 

2.3.1  Method  A 34 

2.3.2  Method  B 35 

2.4  Beamformer  Performance  Measures 35 

2.4.1  Eigenvalues 35 

2.4.2  Sensitivity  Ratio 37 

2.4.3  Directivity  Index 38 

3  Beamforming  Results  40 

3.1     Preliminaries 40 


3.2  Broadside  Results 41 

3.2.1  Method  A 45 

3.2.2  Method  B 52 

3.3  Steered  Results 59 

3.3.1  Method  A 64 

3.3.2  Method  B 64 

4     Summary  and  Conclusions  74 

A   Spectral  Covariance  Matrix  Inversion  78 

B    Dolph-Chebyshev  Weights  81 


List  of  Figures 


2-1     Array  Geometry     14 

2-2     Array  Coordinate  System  and  Dimensions 14 

2-3     Element  Numbering  Scheme 19 

2-4     Spectral  Covariance  Matrix 20 

2-5     System  Architecture  for  a  200  Element  Array 23 

2-6  Projected  Scanning  Vector  Locations.  Every  other  scanning  increment  in 
elevation  and  azimuth  is  shown.  The  shading  illustrates  the  variability  in 
azimuth  of  the  spatial  separation  between  scan  directions  relative  to  the 
spatial  separation  at  broadside  resulting  from  the  uniform  angular  scanning 

implementation 28 

2-7  Null  Placement  Effects.  For  a  source  (solid  line)  incident  at  some  k,  the 
MVDP  beamformer  optimally  places  a  null  in  the  array  response  pattern 

(dashed  line)  in  the  direction  of  the  source  incident  at  k 30 

2-8     Three  Dimensional  Mesh  Plot  of  Penalty  Function 32 

2-9     Overhead  View  of  Penalty  Function     33 

2-10  Two  Dimensional  Cross  Section  of  Penalty  Function 33 

3-1     array  output  for  plane  wave  case 42 

3-2     conventional  beam  pattern  steered  to  0°  elevation  and  90°  azimuth  ....  43 

3-3     conventional  beam  pattern  steered  to  0°  elevation  and  102°  azimuth    ...  44 

3-4     Eigenvalue  Performance  Results  for  Case  11 A 47 

3-5     array  output  for  case  HA,  Amtn  =  10-3     48 

3-6     beam  pattern  for  case  HA,  Amtn  =  10-3 49 


3-7     square  magnitude  of  element  weights,  case  11  A,  Amtn  =  10   "      50 

3-8     corrected  phase  of  element  weights,  case  11  A,  Amtn  =  10~ 51 

3-9     Eigenvalue  Performance  Results  for  Case  9B 53 

3-10  array  output,  case  9B,  Am$n  =  10_ 55 

3-11  beam  pattern,  case  9B,  Amtn  =  10-4 56 

3-12  square  magnitude  of  element  weights,  case  9B,  Amin  —  10-4 57 

3-13  corrected  phase  of  element  weights,  case  9B,  Amin  =  10         58 

3-14  array  output,  case  5B,  Amtn  =  10-4 60 

3-15  beam  pattern,  case  5B,  Amtn  =  10_      61 

3-16  square  magnitude  of  element  weights,  case  5B,  Amin  =  10        62 

3-17  corrected  phase  of  element  weights,  case  5B,  Amtn  =  10-4 63 

3-18  array  output,  steered  method  A  penalty  function 65 

3-19  beam  pattern,  steered  method  A  penalty  function 66 

3-20  square  magnitude  of  element  weights,  steered  method  A  penalty  function  67 

3-21  corrected  phase  of  element  weights,  steered  method  A  penalty  function  .  .  68 

3-22  array  output,  steered  method  B  penalty  function 70 

3-23  beam  pattern,  steered  method  B  penalty  function 71 

3-24  square  magnitude  of  element  weights,  steered  method  B  penalty  function  72 

3-25  corrected  phase  of  element  weights,  steered  method  B  penalty  function  .  .  73 

4-1  Plot  of  First  Sidelobe  versus  Center  Region  Width,  Method  A.  Conven- 
tional main  lobe  widths  in  terms  of  angular  separation  between  first  nulls 

are  approximately  24°  in  azimuth  and  50°  in  elevation 75 

4-2  Plot  of  First  Sidelobe  versus  Center  Region  Width,  Method  B.  Conventional 
main  lobe  widths  in  terms  of  angular  separation  between  first  nulls  are 

approximately  24°  in  azimuth  and  50°  in  elevation 76 

B-l    Dolph-Chebyshev  square  magnitude  element  weights 82 


List  of  Tables 


2.1     Summary  of  Computational  Requirements  for  Different  Scan  Patterns    .  .  27 

3.1  Penalty  Function  Parameters 45 

3.2  Minimum  Eigenvalues  and  Associated  Dynamic  Ranges  Tested  for  Method 

A  Penalty  Functions 46 

3.3  Method  A  Case  Study  Summary 47 

3.4  Minimum  Eigenvalues  and  Associated  Dynamic  Ranges  Tested  for  Method 

B  Penalty  Functions 53 

3.5  Method  B  Case  Study  Summary 54 

A.l    Computational  Times  for  Matrix  Solution  Methods 79 


Chapter  1 


Introduction 


1.1      Background 

The  problem  of  beamforming  with  multidimensional  arrays  is  a  very  important  prob- 
lem with  numerous  applications,  yet  it  is  not  well  understood.  Whereas  one  dimensional, 
equally  spaced  arrays  have  many  formulations  and  solutions  based  on  extensive  temporal 
analysis,  arrays  with  added  dimensions  and  nonuniform  spacing  introduce  unique  prob- 
lems and  possibilities  which  have  no  direct  temporal  equivalent. [l]  In  terms  of  the  bearing 
estimation  problem  which  is  addressed  in  this  investigation  conventional  beamforming 
can  be  described  as  the  shading  or  windowing  of  array  sensor  data  to  obtain  a  desirable 
directional  scan  pattern,  or  beam  pattern.  One  dimensional  windows  are  numerous  and 
well  documented.  [2]  In  two  dimensional  and  higher  dimensional  cases,  fewer  windows  are 
available. [3]  The  majority  of  the  multidimensional  windows  are  based  either  on  a  perfect 
circular  based  symmetric  geometry  (a  perfect  sphere  or  cylinder)  or  on  the  extension  of 
one  dimensional  windows  to  higher  dimensions  which  routinely  requires  uniformly  sam- 
pled apertures. [4,5]  Unfortunately,  all  multidimensional  arrays  are  not  completely  circular 
and/or  do  not  have  regular  sensor  spacing.  The  impact  of  array  geometry,  particularly  in 
three  dimensions,  on  the  spectral  estimation  and  bearing  estimation  problems  is  not  well 
understood.  [6]  Therefore,  the  purpose  of  this  investigation  is  to  develop  a  beam  pattern 
design  procedure  that  can  be  used  on  general  multidimensional  irregular  arrays  which 
includes  the  effects  of  the  array  geometry.    The  design  procedure  uses  a  data-adaptive, 


minimum  variance  method  that  can  be  applied  to  arrays  with  nonuniform  spacing.  The 
minimum  variance  procedure  is  Capon's  High  Resolution  Frequency  Wavenumber  Esti- 
mation Method  more  commonly  known  as  the  Maximum  Likelihood  Method  (MLM).[7] 
The  test  array,  which  has  proposed  applications  in  an  undersea  environment,  consists  of 
a  three  dimensional  sonar  array  located  on  an  arc  of  a  cylindrical  surface.  The  primary 
objective  is  to  design  beam  patterns  by  processing  array  data  with  the  MLM  beamformer. 
The  source  of  array  data  is  the  design  tool  and  is  contained  in  a  special  spectral  covariance 
matrix  termed  a  "penalty  function."  Generation  of  the  penalty  function  is  the  important 
issue  of  the  design  procedure  and  is  the  focus  of  this  investigation. 

1.2      Objectives 

As  previously  stated,  the  primary  objective  is  to  develop  a  beam  pattern  design  pro- 
cedure using  an  MLM  beamformer  that  can  be  used  on  irregular  arrays  where  traditional 
aperture  shadings  are  not  applicable.  This  primary  objective  results  in  two  issues: 

•  can  the  MLM  beamformer  be  effectively  used  in  a  design  role? 

•  can  the  MLM  beamformer  control  sidelobes? 

Both  of  these  issues  are  examined  using  the  penalty  function  concept.  Selection  of  the 
proper  penalty  function  results  in  a  beam  pattern  with  low  sidelobes  when  compared 
to  a  conventional  beam  pattern  generated  from  a  uniform  aperture  shading.  Different 
penalty  functions  result  in  beam  patterns  with  different  sidelobes.  The  selection  of  dif- 
ferent penalty  functions  gives  the  designer  control  over  the  procedure.  In  applying  this 
procedure  the  designer  obtains  a  qualitative  perception  of  the  array's  performance  in 
terms  of  the  array  geometry  and  element  spacing  which  aids  in  the  design  process. 

To  aid  in  the  analysis  and  the  evaluation  of  the  primary  objective,  quantitative  mea- 
sures of  array  performance  in  terms  familiar  to  conventional  beamforming  are  included. 
The  two  performance  measures  relied  upon  are  directivity  index  and  sensitivity  ratio. 
Sensitivity  ratio  is  a  relative  measure  of  white  noise  array  gain  or  array  sensitivity.  The 
quantitative  measure  used  to  describe  the  penalty  function  in  terms  of  array  performance 
is  the  minimum  eigenvalue  of  the  penalty  function.  One  additional  aid  is  a  comparison  of 
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the  shading  obtained  by  the  design  procedure  with  the  aperture  shading  resulting  from  a 
Dolph-Chebyshev  design  for  a  uniform  rectangular  array.  The  quantitative  performance 
measures  and  Dolph-Chebyshev  comparison  study  enable  an  analysis  and  evaluation  based 
on  proven  conventional  methods  and  parameters. 

The  array  application  used  in  this  investigation  deals  with  a  large  array  composed  of 
200  elements  in  a  three  dimensional  sound  field.  Arrays  with  a  large  number  of  elements 
(as  is  the  case  here)  lead  to  computationally  complex  and  time  intensive  array  processing 
requirements;  therefore,  a  significant  issue  concerns  obtaining  an  implementation  which 
reduces  or  limits  the  computational  complexity  and  processing  times.  In  order  to  limit 
the  processing  requirements,  the  aspects  of  directional  (baffled)  hydrophones  or  sensors, 
array  structural  shading  and  the  estimation  process  of  the  spectral  covariance  matrix  are 
not  considered.  The  array  structure  is  transparent  to  sound  from  any  direction  and  the 
spectral  covariance  matrix  is  assumed  to  be  available.  The  procedure  is  based  on  a  plane 
wave  development.  Since  the  MLM  requires  a  matrix  inverse,  even  the  use  of  these  sim- 
plifications still  leads  to  a  time  intensive  process.  In  an  effort  to  further  reduce  processing 
times,  the  investigation  highlights  a  comparison  of  two  matrix  solution  methods: 

•  Gaussian  elimination 

•  a  Toeplitz  bordering  approach  based  on  Levinson's  method.  [8] 

1.3      Overview 

The  remainder  of  this  paper  is  divided  into  three  chapters  and  two  appendices. 
Chapter  2  describes  the  problem  and  is  broken  into  four  additional  sections. 

•  Section  2.1  describes  the  array  geometry. 

•  Section  2.2  defines  the  signal  model,  develops  the  spectral  covariance  matrix  and 
its  properties,  reviews  conventional  beamforming  and  MLM  beamforming,  defines 
array  gain  and  develops  the  spatial  scanning  problem. 

•  Section  2.3  develops  the  penalty  function. 
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•  Section  2.4  develops  and  defines  the  beamformer  performance  measures. 

Chapter  3  contains  the  beamforming  design  results.  The  design  results  are  split  into 
three  sections. 

•  Section  3.1  includes  beamformer  verification  results  which  are  useful  as  points  of 
comparison. 

•  Section  3.2  contains  extensive  results  on  broadside  beam  pattern  designs. 

•  Section  3.3  contains  results  of  steered  beam  pattern  designs. 

Chapter  4  contains  a  summary. 

Appendix  A  contains  timing  data  resulting  from  the  investigation  of  the  Gaussian 
elimination  and  Toeplitz  bordering  matrix  solution  methods. 

Finally,  appendix  B  contains  a  comparison  of  the  aperture  shading  design  results  with 
the  traditional  Dolph-Chebyshev  design  procedure. 
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Chapter  2 


Problem  Development 


2.1      Array  Geometry 

The  array  used  in  this  investigation  has  proposed  applications  on  submersible  bodies. 
It  is  located  on  a  60  degree  arc  of  a  cylindrical  surface.  The  array  is  shown  in  figure  2-1 
and  is  to  be  operationally  oriented  as  depicted.  The  array  has  intrinsic  structure  but  is 
neither  linear  nor  circular.  The  array  is  to  be  applied  in  a  manner  which  reduces  the 
problem  to  processing  only  the  half  space  on  the  same  side  of  the  cylinder  on  which  the 
array  is  located. 

The  array  consists  of  200  elements  arranged  in  ten  rows  (or  lines)  of  twenty  elements 
each.  Considered  separately,  each  line  of  twenty  elements  is  identical.  The  array  is 
symmetric  about  the  x  and  z  axes  for  the  coordinate  system  defined  in  figure  2-2.  This 
structure  leads  to  some  important  results  when  considering  the  spectral  covariance  matrix. 

The  design  wavelength  is  defined  as  the  wavelength  (A)  that  determines  the  element 
spacing.  All  of  the  work  performed  in  the  investigation  of  this  array  is  accomplished  at  the 
design  wavelength;  therefore,  design  wavelength  and  A  are  interchangeable.  The  element 
spacing  is  specified  to  be  one  quarter  of  the  design  wavelength  or  A/4.  All  dimensions  of 
the  array  (including  the  cylinder  radius)  are  related  in  terms  of  the  design  wavelength. 

The  structure  of  the  array  can  be  compared  to  a  simpler  two  dimensional  planar  array. 
This  comparison  is  valuable  because  uniform  planar  arrays  are  well  documented.  In  order 
to  prove  this  comparison  is  justified  to  first  order,  the  extent  the  actual  three  dimonsional 
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Figure  2-1:  Array  Geometry 


0.25A 


0.223A 


#  =  2.15A 


8  =  0.228A 

Figure  2-2:  Array  Coordinate  System  and  Dimensions 
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array  deviates  from  a  true  uniform  planar  array  must  be  determined.  If  the  array  is 
considered  to  be  a  structure  of  ten  20  element  line  arrays  ,  then  the  length  (L)  of  the  line 
arrays  is  4.75A.  If  the  angular  separation  of  the  line  arrays  is  uniform  and  the  arc  length 
of  this  angular  separation  is  approximately  linear  and  equal  to  A/4,  then  the  relationship 
between  R  and  A  can  be  determined  by  equating  the  arc  length  of  the  60  degree  sector  to 
the  arc  length  enscribed  by  the  10  line  arrays.  This  results  in 

-R  =  9-  (2.1) 

3  4  y      ' 


or 


27 
R  =  —  A  (2.2) 

By  applying  some  simple  trigonometry,  the  planar  height  (H)  of  the  array  as  projected 
onto  the  x-z  plane  is  equal  to  the  radius  R,  or  2. 15 A.  The  maximum  deflection  of  the  array 
(or  a  measure  of  the  curvature  from  the  x-z  plane)  is  approximately  equal  to  0.288A.  The 
projected  vertical  spacing  in  the  array  center  is  equal  to  0.25A  and  is  approximately  equal 
to  0.223A  at  the  array  top  or  bottom.  A  pictorial  explanation  of  these  dimensions  is 
illustrated  in  figure  2-2.  The  projected  array  has  uniform  spacing  along  the  axial  (or  z) 
axis  and  nonuniform  spacing  in  the  x  direction  (although  the  spacing  is  nearly  uniform.) 
Therefore,  the  three  dimensional  array  can  be  approximated  by  a  uniform  rectangular 
array  if  the  effects  of  curvature  are  neglected.  The  rectangular  array  approximation 
is  better  in  the  axial  direction;  therefore,  the  expectation  is  that  design  results  in  the 
axial  direction  will  more  closely  resemble  results  of  the  uniform  rectangular  array.  The 
rectangular  array  comparison  provides  a  reasonable  model  which  is  useful  in  establishing 
an  intuitive  feel  for  the  expected  performance  of  the  test  array. 

2.2      Array  Processing  Approach 
2.2.1      Signal  Model 

One  purpose  of  an  array  and  its  associated  processor  is  to  provide  an  estimate  of 
the  directional  spectrum  of  a  space/time  stochastic  process.    The  directional  spectrum 
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information  is  contained  in  the  frequency  wavenumber  function.  To  understand  the  ar- 
ray processing  system  the  signal  model  must  first  be  defined.  The  space/time  stochastic 
processes  being  considered  are  assumed  to  be  temporally  wide  sense  stationary  (defined 
in  reference  [9])  and  spatially  homogeneous  (defined  in  reference  [l].)  The  wide  sense 
stationarity  assumption  implies  that  disjoint  frequency  bands  are  uncorrelated.  The  spa- 
tial homogeneity  assumption  implies  that  disjoint  wavenumber  bands  or  directions  are 
uncorrelated,  i.  e.  no  correlated  multipaths.  Both  the  wide  sense  stationarity  and  spatial 
homogeneity  assumptions  are  not  universally  valid  in  the  ocean  environment,  but  can  be 
closely  approximated  in  practice  by  proper  choices  of  time  and  space  scales.  Although 
time  domain  analysis  is  not  addressed  in  this  study  and  the  actual  physical  dimensions  of 
the  array  are  not  known,  the  wide  sense  stationarity  and  spatial  homogeneity  assumptions 
are  presumed  to  be  valid  based  on  anticipated  time  and  spatial  scales.  The  space/time 
stochastic  processes  are  represented  in  a  Stieltjes  integral  form  [10]  as 

/oo      roo     roo 
/       /       e^"-k')<fX(/,k)  (2.3) 

-oo  J  —oo  J  —oo 

where  t  is  the  time  argument,  z  is  the  vector  representing  the  spatial  position  and  k  is 
the  three  dimensional  wavenumber  vector.1  For  the  wide  sense  stationary  and  spatially 
homogeneous  processes  considered,  X(f,  k)  defines  a  Fourier  transform  representation  for 
a  plane  wave  of  frequency  f  and  wavenumber  vector  k.  Thus  x(t,  z)  is  composed  of  a 
superposition  of  plane  waves  X(/,k).[l]  For  the  three  dimensional  ocean  environment,  / 
and  k  must  satisfy  the  dispersion  relation 

|k|  =  2ir//c  (2.4) 

Equation  2.4  results  from  the  three  dimensional  wave  equation  where  c  is  the  speed  of 
sound  in  water.  This  dispersion  relation  corresponds  to  free  space  propagation.  The 
space/time  covariance  function  K(t,  t—T,  z,  z— Az)  is  the  cross  correlation  of  the  stochastic 
process  at  times  t  and  t  —  r  and  spatial  positions  z  and  z  —  Az.  All  processes  are  assumed 


1  Lower  case  letters  are  used  to  depict  temporal  representations  while  upper  case  letters  depict  transform 
representations.  Similarly,  upper  case  bold  face  letters  represent  matrices  and  lower  case  bold  face  letters 
represent  vectors. 
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to  be  zero  mean  stochastic  processes.  Kz(t,  t  -  t,  z,  z  -  Az)  is  determined  by 

Kx(t,t-r,z,z-  Az)  =  E[x(t,z)x*(t-r,z-  Az)]  (2.5) 

where  *  denotes  complex  conjugate  and  E  represents  the  ensemble  expectation  operator. 
The  frequency  domain  representation  of  the  space/time  covariance  function  is  the  spectral 
covariance  function  Sx(f,  z,z  —  Az)  which  is  defined  as 

/oo 
Kx(t,  t-T,z,z-  Az)e~ j2*fTdr  (2.6) 

-oo 

Sz(f,  z,z  —  Az)  represents  the  cross  spectra  between  positions  (corresponding  to  array 
elements)  z  and  z  —  Az.  The  frequency  wavenumber  representation  of  the  space/time 
stochastic  process,  called  the  frequency  wavenumber  function  P(f,  k),  is  defined  as 

P(f,  k)  =  J  J  |~  8.{ft  z,  z  -  AzykAzd(Az)  (2.7) 

Since  /  and  |k|  are  related  by  equation  2.4,  fixing  /  reduces  the  problem  to  determining 
the  direction  of  propagation  of  |k|.  This  is  the  approach  used  throughout  this  study;  / 
is  fixed  and  then  the  half  space  in  which  the  array  resides  is  scanned  to  determine  the 
direction  of  propagation  parameterized  in  terms  of  an  bearing  in  azimuth  and  an  elevation 
angle.  The  estimate  of  the  direction  of  propagation  or  directional  spectrum,  embedded  in 
the  function  P(f,  k),  is  the  desired  output  of  the  array  processor. 

2.2.2      Spectral  Covariance  Matrix 

An  examination  of  the  spectral  covariance  matrix  is  included  at  this  point  because  it 
aids  in  the  understanding  of  the  beamforming  procedure  and  also  reveals  some  pertinent 
characteristics  of  the  matrix  operations  which  are  encountered. 

The  spectral  covariance  matrix  can  be  decomposed  into  two  components,  one  coherent 
and  one  incoherent,  by  setting  (dropping  arguments  for  brevity)2 

Sx  =  a|ddt  +  a2NQ  (2.8) 


2*  represents  complex  conjugate  transpose 
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where 


,yk«Z] 


JK-*n 


(2.9) 


is  the  geometric  phase  or  direction  vector  of  a  plane  wave  signal  propagating  with  wavenum- 
ber  vector  ka  and  signal  strength  erf.  Q  represents  the  normalized  (trace  Q  =N)  noise 
cross  spectral  density  matrix  of  element  outputs  with  noise  strength  a2N.  If  the  noise 
consists  solely  of  uncorrelated  components,  then  Q  reduces  to  the  identity  matrix  I. 

The  coherent  input  places  the  restrictions  on  /  and  kB  found  in  equation  2.4.  The 
coherent  input  used  in  this  study  consists  of  plane  waves  which  are  modeled  at  time  t  and 
position  z  as 

x(t,z)=a(/yW-k'z)  (2.10) 

where  a(f)  is  a  frequency  dependent  amplitude  function.  The  coherent  component  of  the 
spectral  covariance  matrix  resulting  from  the  plane  wave  input  (determined  by  applying 
equations  2.5  and  2.6)  is 

a2sdS  =  <rf  e->'k'A^       V  ij  (2.11) 

where  cr|  is  frequency  dependent  and  AztJ-  =  zt-  —  z»-  corresponds  to  the  vector  difference 
between  the  t  and  j  elements.  The  exponential  argument  ks  •  Azt)J-  represents  the 
phase  delay  due  to  propagation  of  the  plane  wave  between  the  i     and  j     elements. 

The  first  observation  deals  with  the  structure  of  a  matrix  formed  by  equation  2.8. 
Matrices  of  this  form  are  Hermitian  symmetric  or  simply  Hermitian.  A  Hermitian  matrix 
is  defined  by  S  =  S^.  The  second  observation  results  after  examining  equation  2.11. 
The  coherent  input  is  dependent  only  on  the  difference  in  element  positions.  This  spatial 
homogeneity  leads  to  an  embedded  Toeplitz  structure.  A  Toeplitz  matrix  is  defined  as  a 
matrix  whose  ijth  element  s(i,  j)  is  a  function  of  (i-j)  and  thus  has  identical  elements  along 
the  main  diagonal  and  each  of  the  subdiagonals.  (Reference  [11]  contains  an  excellent 
summary  of  Hermitian  and  Toeplitz  matrix  properties.)  To  fully  see  the  structure  of 
the  matrix  formed  by  equation  2.11  for  the  test  array,  consider  an  element  identification 
scheme  where  the  elements  are  numbered  consecutively  along  each  successive  row  (or  line) 
starting  at  the  array  top.    This  numbering  scheme  is  shown  in  figure  2-3.   The  spectral 
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61   X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  80 

81  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  X  100 
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161  XXXXXXXXXXXXXXXXXXXX  i8o 

181  xxxxxxxxxxxxxxxxxxxx  200 


Figure  2-3:  Element  Numbering  Scheme 

covariance  matrix  resulting  from  this  numbering  scheme  is  shown  in  figure  2-4  where  the 
matrix  is  partitioned  to  indicate  specific  sub-matrices.  Each  sjfcj  corresponds  to  the  cross 
spectral  covariance  of  the  kth  element  with  the  Ith  element.  Each  sub-matrix  corresponds 
to  the  spectral  covariance  formed  between  elements  located  in  two  of  the  line  arrays.  The 
main  diagonal  of  sub-matrices  is  the  result  of  forming  the  spectral  covariance  of  each  line 
array  with  itself;  the  main  diagonal  sub-matrices  are  all  identical  because  each  line  array 
is  identical  and  the  coherent  input  is  perfectly  correlated. 

Sx  can  also  be  represented  as  a  two  level  matrix  where  to  the  first  order  level 

Si,i      Si(2     ••*      Si_io 

S2,l       S2,2 

(2.12) 

Sio.i  Sio.io 

Each  sub-matrix  St)J-  is  the  cross  covariance  between  the  ith  and  jth  line  arrays  and 
is  defined  in  figure  2-4  by  each  partition.  For  example  S12  is  the  sub-matrix  in  the 
partitioned  first  row  and  second  column  of  figure  2-4  and  is  formed  from  the  line  array 
with  elements  1  through  20  and  the  line  array  with  elements  21  through  40.  The  spectral 
covariance  matrix  has  several  distinctive  characteristics. 

•  Sx  is  NxN  (or  200x200  for  the  test  array.) 

•  Each  sub-matrix,  S,y,  is  Toeplitz  and  Hermitian 

•  Sjb,  to  the  first  order  level,  is  Hermitian  but  not  Toeplitz  because  of  the  vertical 
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Figure  2-4:  Spectral  Covariance  Matrix 
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structure  of  the  array.  However,  the  main  diagonal  is  a  constant  (to  the  first  order 
level.) 

•  The  addition  of  the  incoherent  term  will  not  change  the  embedded  Toeplitz  structure 
of  the  sub-matrices. 

Matrices  with  the  form  of  Sx  are  frequently  classified  as  a  two  level  matrices.  Ref- 
erences [12]  and  [13]  include  a  complete  discussion  on  the  classification  of  matrices  in 
this  manner.  Sx  in  this  case  is  classified  as  a  Hermitian-Toeplitz  two  level  matrix,  but 
by  permuting  corresponding  rows  and  columns,  the  matrix  can  be  transformed  into  a 
Toeplitz-Hermitian  two  level  matrix  where  the  first  level  is  Toeplitz  and  the  second  level 
sub- matrices  are  Hermitian.  This  is  important  because  there  currently  exist  efficient  and 
documented  algorithms  which  can  invert  a  Toeplitz-Hermitian  two  level  matrix. [13]  As  is 
shown  later  the  spectral  covariance  matrix  needs  to  be  inverted.  Details  of  the  matrix 
inversion  selection  process  are  included  in  Appendix  A.  (This  selection  process  played  a 
very  important  role  in  this  investigation;  the  time  intensive  inversion  process  had  a  direct 
relationship  to  the  processing  time.) 

2.2.3     MLM  Development 

In  the  context  of  this  investigation  array  processing  and  beamforming  are  synony- 
mous. There  is  much  literature  on  beamforming,  both  conventional  and  data  adaptive. 
Complete  reviews  of  processing  methods  are  readily  available.  [14]  Since  the  purpose  of  this 
study  is  to  develop  a  beam  pattern  design  procedure  that  can  be  applied  to  nonuniform 
multidimensional  arrays,  the  choice  of  usable  methods  is  limited  to  the  MLM  beamformer. 
The  procedure  can  be  more  accurately  described  as  a  minimum  variance,  distortionless 
processor  (sometimes  referred  to  as  MVDP).[6]  From  this  point  on,  the  beamformer  will 
be  referred  to  as  the  MVDP  because  it  more  accurately  describes  the  processor.  What 
follows  is  a  development  of  conventional  beamforming  and  the  MVDP. 

w(/|kr)  is  the  aperture  weighting  function  with  the  processor  steered  to  k^.3  For  the 
sampled  aperture  consisting  of  an  N  element  array,  the  »"*"  component  of  w  is  the  complex 


'The  development  is  the  same  as  used  in  reference  [l]. 
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conjugate  of  the  ith  element's  weighting  (ty*). 


w(/  |  kr)  - 


<(/|kr) 


(2.13) 


The  aperture  weighting  function,  or  weighting  vector,  can  also  be  considered  to  be  an  n- 
dimensional  taper  or  shading  function  for  an  n-dimensional  array.  It  should  be  observed 
that  w  is  dependent  on  the  steer  direction  and  can  be  complex.  The  scanning  vector 
e(kr)  steered  to  ky  is  defined  as 


e(kr) 


jlLT-%N 


(2.14) 


where  again  zt-  is  the  i  element's  position  vector  of  an  N  element  array.  The  scanning 
vector  can  be  interpreted  as  the  plane  wave  phase  delays  which  must  be  applied  to  the 
array  elements  to  steer  the  array  in  the  direction  of  ky.  The  scanning  vector  also  has  an 
interpretation  of  being  a  replica  field  or  the  field  at  the  array  produced  by  a  source  at 
some  position  in  the  water  column. 

The  output  of  the  array  beam  steered  to  ky  is  defined  in  conventional  beamforming 


as 


N 


Y(/,kr)  =  w*X=  X>(/|kT)X(/,zO 


(2.15) 


t=i 


This  system  architecture  is  depicted  in  figure  2-5.  The  array  response  pattern  W(f,  k|kj) 
is  given  by 


N 


W,k|kT)  =  w*e  -  I>(/Ifcryk'** 


(2.16) 


t=i 


and  represents  the  response  of  the  array  (steered  in  the  direction  of  ky)  over  the  aperture 
to  a  plane  wave  signal.  If  k  is  fixed  and  the  function  W(f,  k|kj)  is  evaluated  over  a  spatial 
region,  then  the  beam  pattern  B  is  obtained  where 

B(/,Mkr)  =  ^(/,k|kT)|k=2,a(^}  (2.17) 

where  A  is  the  wavelength  and  a(9,<f>)  is  a  unit  vector  specified  in  a  spherical  geometry. 
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Figure  2-5:  System  Architecture  for  a  200  Element  Array 

Because  P(/,k)  is  a  measure  of  power,  this  leads  to  a  representation  of  the  response  of 
the  array  as  a  function  of  two  array  sensor  locations.  This  will  result  in  a  second  moment 
representation  for  the  array  output  and  requires  a  second  moment  input.  The  array  input 
is  the  spectral  covariance  function  Sx(f,  z,,  Zy)  as  defined  by  equation  2.6.  The  output  of 
the  beam  steered  to  kj  is  the  power  density  spectrum  jSy(/|kx)  and  is  given  by 


5v(/|kr)=wT(/|kT)SI(/,z,-,zy)w(/|kr) 


(2.18) 


iSy(/|kr),  as  defined  by  equation  2.18,  is  the  output  of  a  conventional  beamformer  and 
is  the  estimate  of  P(f,  k)  with  the  array  steered  to  ky.  The  MVDP  is  obtained  by 
minimizing  the  response  of  equation  2.18  (minimum  variance)  subject  to  unity  gain  in  the 
look  direction  (distortionless  response.)  The  MVDP  solutions  are  [1,7,15] 

S;1e(kr) 


WMVDp{f\kT) 


et(kr)Sx1e(kT) 


V)  rf  Vile  ^_  et(k)s;1e(kr) 

et(kT)SI1e(kT) 

eT(kT)Sx1e(kT) 


(2.19) 

(2.20) 
(2.21) 
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The  MVDP  method  adaptively  finds  different  aperture  weightings  for  different  input  signal 
fields  such  that  sidelobes  are  adjusted  to  reject  interfering  coherent  noise  sources  and 
the  incoherent  noise.  This  method  is  optimum  in  the  sense  that  it  minimizes  the  array 
response  over  the  aperture  of  the  noise  field  (the  input  signal  not  in  the  beam)  by  design 
and  application  of  a  spatial  filter.  This  spatial  filter  estimates  the  power  in  a  plane  wave 
propagating  with  frequency  /  and  wavenumber  vector  kr  in  the  presence  of  interfering 
plane  waves  and  sensor  noise. [16]  To  fully  define  the  directional  spectrum,  the  signal  field 
must  be  scanned  over  the  appropriate  angular  space  (which,  for  this  study,  is  the  half  space 
in  which  the  test  array  is  located.)  If  only  uncorrected  noise  is  present,  then  uniform 
weighting  with  the  conventional  beamformer  is  optimal. 

2.2.4     Array  Gain 

Array  gain  is  the  improvement  in  signal  to  noise  ratio  due  to  beamforming.  This  is 
a  result  of  finite  beam  widths  rejecting  noise  and  signals  not  in  the  beam.  Ideally,  for  the 
signal  of  interest,  signal  rejection  is  small  while  the  noise  rejection  is  large  resulting  in  a 
large  array  gain.  The  array  gain  as  such  is  defined  as 

„       lwtd!2  ,        * 

G  =  X—r-±-  2.22 

wTQw  v         ' 

which  is  simply  the  ratio  of  the  signal  response  to  the  noise  response.  In  the  context  of 
the  MVDP  processor,  the  optimum  array  gain  is  [1,17] 

GoVt  =  etQ^e  (2.23) 

For  the  case  of  spatially  white  or  sensor  noise  with  omnidirectional  sensors,  the  array  gain 
becomes  the  white  noise  gain  or 

Gw  =  ^—^-  <  N  (2.24) 

w'w  v         ' 

where  N  is  the  number  of  sensors.  Equality  occurs  when  the  processor  is  perfectly  matched 
and  there  is  uniform  weighting.  The  processor  is  defined  to  be  perfectly  matched  to  the 
signal  directional  characteristics  when  e  =  d  where  d  is  defined  by  equation  2.9. [17] 
Mismatch  occurs  whenever  e  /  d  and  may,  for  example,  be  caused  by  non-plane  wave 
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propagation,  errors  associated  with  the  sensors,  sampling,  quantization,  and  steering  off 
of  the  signal  direction.  When  the  errors  due  to  mismatch  are  uncorrelated,  the  sensitivity 
(Sw)  of  the  array  gain  due  to  signal  mismatch  is  related  to  Gw  by  [18] 

Sw=-±-  (2.25) 

Used  in  this  context,  Sw  or  Gw  is  a  measure  of  robustness.  Using  v?mvdp  calculated 
in  equation  2.19  for  a  perfectly  matched  processor  with  Sz  =  I,  Gw  equals  N  and  the 
beamformer  is  robust.  (yfMVDP  f°r  this  case  is  a  normalized,  uniform  weighting.)  Any 
deviation  from  N  indicates  a  degradation  in  beamformer  robustness. 

One  additional  expression  for  the  array  gain  which  turns  out  to  be  useful  relates  G  in 
terms  of  eigenvalues  and  eigenvectors  of  Q.  The  expression  is 

G  =  E«  (2.26) 

where  A,-  and  vt-  are  the  eigenvalues  and  eigenvectors,  respectively.  [18]  This  expression 
shows  that  G  is  dominated  by  small  ratios  of  |d'v,-|2/Af.  This  effect  is  not  well  understood 
because  when  At-  is  small,  d'V|  is  usually  small  also. 

2.2.5     Spatial  Scanning 

The  remaining  feature  to  be  specified  is  the  spatial  scanning  problem.  The  spa- 
tial field,  as  previously  defined,  is  the  half  space  in  which  the  array  is  located.  To  fully 
define  the  directional  spectrum,  the  spatial  field  must  be  scanned  or  sampled.  Equa- 
tion 2.14  defines  the  method  yet  does  not  answer  the  question,  "how  finely  must  this  field 
be  sampled?"  The  spatial  sampling  interval  is  defined  in  terms  of  Nyquist's  sampling 
requirements.  This  requires  sampling  at  one  half  the  natural  beam  width  of  the  array 
which  is  equivalent  to  sampling  at  twice  the  highest  frequency  of  the  field  or  spatially  at 
A/2.  Three  factors  must  be  considered  when  determining  the  sampling  requirements. 

1.  As  described  in  reference  [19],  finer  sampling  results  in  a  visual  display  with  signifi- 
cantly more  content.  This  additional  content  is  valuable  when  evaluating  hard  copy 
or  video  outputs  of  the  directional  spectrum. 
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2.  If  the  natural  beam  width  is  determined  by  assuming  the  array  is  approximately 
rectangular,  then  error  will  be  introduced  into  the  values  obtained. 

3.  The  MVDP  beamformer  is  susceptible  to  serious  signal  suppression  effects  due  to 
signals  arriving  form  directions  between  the  processor  beams.  Minimizing  this  effect 
requires  more  closely  spaced  scan  directions.  [17] 

Consideration  of  these  three  factors  leads  to  the  selection  of  a  scanning  increment  of 
1/4  the  natural  beam  width  (BW).  Verifying  this  is  the  proper  decision  requires  the 
evaluation  of  1/2  and  1/8  BW  scanning  increments.  There  is  one  additional  consideration 
when  determining  the  proper  scanning  increment;  it  involves  the  computational  burden. 
To  illustrate  this,  the  natural  beam  width  of  the  array  must  be  determined. 

For  the  test  array,  BW  is  resolved  into  two  directions  called  elevation  (also  known 
as  depression/elevation  or  D/E)  and  azimuth.  These  directions  arise  naturally  when 
anticipating  an  operational  scenario;  elevation  is  an  angular  measure  off  of  a  horizontal 
plane  (+90  degrees  to  -90  degrees)  and  azimuth  is  an  angular  measure  fore  and  aft  of  a 
vertical  plane  through  the  center  of  the  array  (0  degrees  forward  to  180  degrees  aft).  To 
determine  BW  the  assumption  that  the  array  is  approximately  uniform  and  rectangular 
is  applied.  BW  for  a  uniform  rectangular  array  in  the  x  direction  is  [20] 

BWZ  =  A  (2.27) 

Lit 

For  the  azimuthal  direction,  Lazimuth  =  4.75A  and 

BWazimuth  =  0.067*  m  12°  (2.28) 

This  leads  to  an  azimuthal  scanning  increment  of  3  degrees.  In  elevation  (where  the 
rectangular  array  approximation  is  least  valid),  Leievation  —  2. 15 A  and 

2° 
BWeUvation  =  0.148*  «  26-  (2.29) 

This  leads  to  a  scanning  increment  in  elevation  of  6|  degrees. 

The  scanning  increment  is  implemented  using  uniform  angular  separation.  The  scan- 
ning system  is  best  described  using  a  global  coordinate  system  of  latitude  and  longitude 
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scan  pattern 

Le 

La 

L 

1/8  BW 

54 

118 

6372 

1/4  BW 

27 

59 

1593 

1/2  BW 

14 

30 

420 

Table  2.1:  Summary  of  Computational  Requirements  for  Different  Scan  Patterns 

where  latitude  corresponds  to  azimuth  and  longitude  corresponds  to  elevation.  The  pro- 
cedure implementing  the  spatial  sampling  is  to  choose  La  latitudes  separated  by  the  az- 
imuthal  scanning  increment.  Then  along  each  latitude,  Le  longitudes  are  chosen  separated 
by  the  elevation  scanning  increment.  The  total  number  of  scans,  L,  is 


Lt  —  LieL/a 


(2.30) 


Steers  in  azimuth  (latitude)  to  0°  and  180°  are  independent  of  the  elevation  angle.  To 
eliminate  duplicate  scans  and  place  a  scanning  direction  at  90°  azimuth,  a  3°  offset  is 
incorporated  at  0°  and  180°  azimuth.  This  pattern  results  in  La  —  59.  To  place  a 
scanning  direction  at  0°,  a  Z\  offset  is  incorporated  at  +90°  and  —90°  elevation;  this 
results  in  Le  =  27.  Therefore,  the  total  number  of  scans  is  L  =  1593.  The  scanning 
pattern  as  projected  onto  a  plane  surface  is  shown  in  figure  2-6. 

Now  consider  the  computational  requirements  of  this  scanning  pattern.  To  define 
the  directional  spectrum  equation  2.21  must  be  calculated  L  times,  once  for  each  scan 
direction.  For  an  N  element  array,  equation  2.21  represents  an  N  by  N  system  of  equations. 
Solving  an  N  by  N  system  of  equations  L  times  is  a  time  intensive  process.  Likewise,  to 
determine  the  beam  pattern,  equation  2.20  must  be  calculated  L  times.  Fortunately,  the 
N  by  N  system  of  equations  only  needs  to  be  solved  once  for  the  array  steer  direction. 
Table  2.1  lists  the  scanning  requirements  for  1/8  BW ,  1/4  BW  and  1/2  BW  scanning 
patterns.  Based  on  these  computational  requirements  the  desire  is  to  sample  the  spatial 
field  at  the  fewest  number  of  points  consistent  with  adequate  beamformer  performance. 

The  scanning  pattern  used  is  not  ideal.  The  algorithm  implementing  this  pattern  is 
very  simple;  with  simple  modifications  1/2  BW  and  1/8  BW  scanning  patterns  can  be 
tested.  However,  the  use  of  uniform  angular  scanning  as  implemented  does  not  translate 
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spatial  separation 
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Figure  2-6:  Projected  Scanning  Vector  Locations.  Every  other  scanning  increment  in 
elevation  and  azimuth  is  shown.  The  shading  illustrates  the  variability  in  azimuth  of  the 
spatial  separation  between  scan  directions  relative  to  the  spatial  separation  at  broadside 
resulting  from  the  uniform  angular  scanning  implementation. 
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into  uniform  spatial  scanning.  Whereas  the  spatial  scanning  is  uniform  in  azimuth,  the 
scanning  in  elevation  is  not  spatially  uniform,  particularly  in  regions  close  to  0°  and 
180°  in  azimuth  where  the  sampling  in  elevation  is  very  dense.  Figure  2-6  illustrates 
the  differences  in  spatial  scanning  increment  in  elevation  relative  to  the  broadside  spatial 
scanning  increment.  If  this  method  of  scanning  is  to  be  practically  implemented,  the 
optimum  way  is  to  rotate  the  pattern  90°  and  limit  the  elevation  scanning  region  to  the 
prime  elevation  angles  of  ±60°  where  no  structural  shading  occurs.  This  results  in  fairly 
uniform  spatial  scanning,  but  over  a  limited  spatial  region. 

In  summarizing  the  spatial  scanning  problem,  a  simple  method  that  suits  the  purposes 
of  the  idealized  synthetic  investigation  of  the  test  array  is  presented.  The  scanning  pattern 
occurs  in  uniform  angular  increments  of  3°  with  59  steers  in  azimuth  and  6|  with  27 
steers  in  elevation  for  a  total  of  1593  steers.  The  spatial  increment  corresponds  to  1/4 
BW  scanning.  The  simple  scanning  method  is  easily  modified  to  allow  the  use  of  different 
scanning  increments,  but  does  not  provide  uniform  spatial  sampling. 

2.3      Beam  Pattern  Design  Method 

In  attempting  to  design  beam  patterns,  the  objective  is  to  reduce  the  sidelobe  levels 
compared  to  sidelobe  levels  characteristic  of  traditional  shading  methods.  The  beam 
pattern  design  procedure  is  implemented  by  the  use  of  a  "penalty  function,"  so  named 
because  the  function  penalizes  beam  pattern  lobes  in  regions  of  high  directional  intensity 
while  enhancing  the  lobes  in  regions  of  low  directional  intensity.  To  understand  this 
function,  the  actual  processes  undertaken  by  the  MVDP  beamformer  must  be  understood 
in  a  qualitative  manner. 

The  MVDP  beamformer  designs  optimum  beam  patterns  by  minimizing  the  energy  in 
the  beam  subject  to  the  constraint  of  unity  gain  in  the  scan  direction.  The  beam  energy 
is  represented  as 

S(f\kT)=  J  P(f,\,)\W\2dk  (2.31) 

which  is  interpreted  as  the  integral  over  all  wavenumber  space  of  the  frequency  wavenum- 
ber  function  weighted  by  the  square  magnitude  of  the  array  response  pattern  (or  beam 
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Figure  2-7:  Null  Placement  Effects.  For  a  source  (solid  line)  incident  at  some  k,  the 
MVDP  beamformer  optimally  places  a  null  in  the  array  response  pattern  (dashed  line)  in 
the  direction  of  the  source  incident  at  k. 

pattern).  The  MVDP  beamformer  designs  array  response  patterns  based  on  the  frequency 
wavenumber  function  of  the  incident  field  in  order  to  minimize  the  beam  energy.  If  a  dis- 
crete source  not  propagating  in  the  look  direction  is  incident  (resulting  in  large  P(f,  k) 
for  that  direction),  then  \W\2  is  reduced  to  minimize  the  beam  energy.  The  MVDP  beam- 
former  optimally  places  nulls  in  the  directions  of  off-look  direction  sources. [15]  This  null 
placement  effect  is  illustrated  in  figure  2-7.  If  diffused  directional  sources  are  located  in  a 
finite  sector,  then  a  lower  array  response  occurs  in  the  direction  of  the  sector.  However, 
the  placement  of  deep  nulls  results  in  higher  overall  sidelobe  levels  which  increases  the 
effects  of  sensor  noise,  the  other  source  term  which  contributes  to  total  beam  energy. 
Conversely,  minimizing  the  effects  of  sensor  noise  by  lowering  overall  sidelobe  levels  re- 
duces the  effects  of  nulls  and  increases  the  beam  energy  from  off-look  direction  sources. 
The  effects  due  to  discrete  directional  sources  and  sensor  noise  are  compromised  until  the 
minimum  beam  energy  is  obtained. 

Since  the  desired  effect  is  a  reduction  in  sidelobe  levels,  the  primary  component  of 
the  penalty  function  consists  of  a  field  of  diffuse  directional  sources  of  varying  intensity. 
The  field  of  directional  sources  is  modeled  similarly  to  isotropic  noise  which  is  directional 
noise  incident  uniformly  from  all  directions  with  uniform  intensity.  The  penalty  function 
directional  source  component  consists  of  a  modified  isotropic  noise  field,  incident  uniformly 
from  all  directions  but  with  nonuniform  intensity.  In  terms  of  equation  2.31,  the  penalty 
function  (defined  by  P(f,  k))  in  the  directions  of  beam  pattern  sidelobes  has  high  intensity 
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while  in  the  direction  of  the  beam  pattern  main  lobe,  has  low  intensity  relative  to  the  high 
intensity  region.  A  three  dimensional  mesh  representation  of  the  penalty  function  is  shown 
in  figure  2-8.  In  regions  of  high  intensity,  \W\2  is  reduced.  In  regions  of  low  intensity,  \W\2 
can  be  large  because  the  weighted  integrated  response  is  small.  If  the  depressed  sector  of 
low  intensity  directional  sources  is  centered  about  the  steer  direction  and  encompasses  an 
area  equivalent  to  the  main  beam  region,  the  main  beam  of  the  resulting  array  response 
pattern  will  be  large  (supported  by  the  unity  gain  constraint  in  the  look  direction)  relative 
to  the  sidelobe  region  which  is  suppressed  due  to  the  uniformly  distributed  high  intensity 
sources  located  in  that  region.  The  net  effect  is  the  desired  result:  a  beam  pattern 
with  suppressed  sidelobes.  Figure  2-9  shows  an  overhead  view  of  the  penalty  function 
centered  about  the  steer  direction.  (The  penalty  function  is  always  centered  about  the 
steer  direction.)  As  shown,  the  penalty  function  consists  of  three  regions. 

1.  The  center  region,  defined  by  width  in  elevation  (A0C)  and  width  in  azimuth  (A^c), 
represents  the  low  intensity  region 

2.  The  transition  region,  defined  by  width  in  elevation  (A6t)  and  width  in  azimuth 
(A<£f),  serves  as  the  interface  between  low  and  high  intensity  regions 

3.  The  surrounding  background  region  of  high  intensity 

The  issues  pertinent  to  penalty  function  design  are  the  shapes  and  widths  of  the  center  and 
transition  regions,  the  amount  of  stabilization,  and  effects  of  the  steer  direction.  Several 
of  these  parameters  are  discussed  now;  others  are  addressed  in  a  later  section. 

The  method  used  to  input  the  penalty  function  into  the  array  processor  must  take  the 
form  of  a  spectral  covariance  matrix.  Using  the  same  concepts  as  found  in  equations  2.8 
and  2.11,  the  penalty  function  PF  is  defined  as 

PF  =  a2NI  +  SN  (2.32) 

where  aN  is  the  sensor  noise  strength  serving  as  a  stabilization  factor  and  Spi  is  the 
directional  component  representing  a  modified  isotropic  noise  field  in  terms  of  a  sum 
of  plane  waves.    The  plane  waves  model  this  directional  field  at  one  half  the  canning 
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Figure  2-8:  Three  Dimensional  Mesh  Plot  of  Penalty  Function 
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Figure  2-10:  Two  Dimensional  Cross  Section  of  Penalty  Function 
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increment.  S]y  is  defined  as 

■k       -k  200  200 

S,=  EE  S(/ A  </>)  E  E  e'T^W*  (2.33) 

^=0  0=0  t=lj=l 

where  S(f,9,<f>)  is  a  spatial  weighting  function  acting  on  the  plane  waves  and  a.r(9,(j>)  is 
a  unit  directional  vector  defining  the  direction  of  propagation.  ]C£=o  Y,e=o  represents  the 
sum  over  the  scanning  space  consisting  of  the  half  space  in  which  the  array  is  located.  S^ 
is  simply  the  sum  of  plane  waves  spatially  weighted  to  form  the  distinct  sectors  shown 
in  figure  2-9.  A  two  dimensional  cross  section  of  S(f,6,<f>)  is  shown  in  figure  2-10.  This 
illustrates  the  regions  of  high  and  low  intensity.  Two  methods  implementing  the  penalty 
function  are  now  examined. 

2.3.1     Method  A 

Method  A  uses  the  direct  application  of  null  placement  to  reduce  beam  pattern 
sidelobes.  The  high  intensity  directional  noise  field  is  formed  with  unit  intensity  while 
the  low  intensity  region  is  formed  with  zero  intensity  (or  no  directional  field).  The  sensor 
noise  is  then  used  to  stabilize  the  process.  The  resulting  directional  noise  field  causes 
the  sidelobes  to  be  suppressed  while  the  main  lobe  level  is  maintained  by  the  unity  gain 
constraint. 

Referring  to  figure  2-10,  method  A  is  implemented  by  setting  Smax(f,8,<f>)  =  1  and 
A/>  =  1.  In  terms  of  equation  2.33,  this  requires 

SA,center{f,Oi4>)=0  (2.34) 


{0-0o)\       ,    {4>-4>o? 
-^-}exp[ - 


SA,transition(f,  0, 4>)  =  exp[-  K\\o)   }exp[-  ^-g^   ]  (2.35) 


SA,background{f,6,<f>)  =  1  (2.36) 

This  produces  flat  center  and  background  regions  interfaced  by  a  Gaussian  shaped  tran- 
sition region.  60  and  <f>0  are  the  directions  representing  the  appropriate  center/transition 
region  boundary.  For  a  fixed  number  of  increments  in  the  transition  region,  the  region 
dimensions  and  roll  off's  can  be  controlled  by  adjusting  the  variance  (aj  for  elevation  and 
o\  for  azimuth)  of  the  Gaussian  function. 
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2.3.2      Method  B 

Method  B  uses  the  opposite  of  null  placement  to  reduce  sidelobes.  Instead  of  placing 
nulls  in  the  sidelobe  region,  anti-nulls  are  placed  in  the  main  lobe  region.  This  is  accom- 
plished by  forming  the  high  intensity  region  with  zero  intensity  (no  directional  sources) 
while  forming  the  low  intensity  region  with  negative  intensity  sources.  The  low  intensity 
region  still  has  lower  intensity  relative  to  the  high  intensity  region.  The  effect  of  sidelobe 
suppression  is  indirectly  accomplished  with  method  B.  The  negative  intensity  (or  low  in- 
tensity) region  can  support  very  large  main  beam  levels  since  the  integrated  beam  energy 
from  equation  2.31  is  small.  This  effectively  enhances  the  main  lobe.  The  sidelobes  are 
suppressed  when  the  unity  gain  constraint  scales  the  large  main  lobe.  The  sensor  noise  is 
again  used  to  stabilize  the  process;  however,  greater  stabilization  is  required  for  method 
B  than  is  required  for  method  A. 

Referring  to  figure  2-10,  method  B  is  implemented  by  setting  Smax(f,6,(f>)  =  0  and 
Ap  —  1.  In  terms  of  equation  2.33,  this  requires 

SB.center(f,0,<f>)  =  -l  (2.37) 

SB,tran.iUon{f,  *,  *)  =  -exp[-  £z^£]CTp[-  ^  "  ^  ]  (2.38) 

SB,background{f>QA)  =  °  (2.39) 

The  region  shapes  and  dimensions  are  controlled  in  the  same  manner  as  in  method  A. 

In  order  to  discuss  the  other  issues,  quantitative  measures  which  enable  an  evaluation 
of  those  relevant  issues  must  be  determined.  In  particular,  measures  are  needed  to  evaluate 
the  effects  of  region  widths  and  the  amount  of  stabilization  required. 

2.4      Beamformer  Performance  Measures 

2.4.1     Eigenvalues 

Equation  2.21  requires  that  some  specific  restrictions  be  placed  on  the  spectral  covari- 
ance  matrix  of  the  penalty  function.  The  most  significant  requirement  is  that  Sx  must  be 
invertible.  For  a  meaningful  problem,  this  is  ensured  by  requiring  Sx  be  positive  definite. 
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The  approach  taken  to  ensure  this  requirement  is  satisfied  centers  around  the  eigenvalues 
of  Sx. 

In  section  2.2.2,  it  was  determined  that  Sx  is  Hermitian  which  results  in  some  special 
eigenvalue  properties.  All  eigenvalues  in  this  investigation  are  calculated  using  the  Hermi- 
tian matrix  routines  available  in  the  EISPACK  eigensystem  package. [21]  Since  eigenvalues 
of  Hermitian  systems  are  always  real  [22],  a  simple  characterization  of  Sx  is  obtained.  This 
simple  characterization  has  physical  interpretation  which  is  now  pursued. 

Consider  the  eigenvalue  problem 

Sxv,-  =  A,v,  (2.40) 

where  A,-  is  the  i**  eigenvalue  and  v,  is  the  eigenvector  corresponding  to  A,-.  Sx  can  be 
characterized  as  a  measure  of  the  energy  of  the  received  signals  (similar  to  the  eigenvalue 
interpretation  of  array  gain.)  This  leads  to  an  energy  interpretation  of  eigenvalues  where  A,- 
is  a  measure  of  the  energy  projected  in  the  direction  of  vt.  Small  A,-  represent  eigenvectors 
with  low  intensity.  Low  energy  levels  are  associated  with  eigenvectors  corresponding  to  the 
low  intensity  region  of  the  penalty  function.  Again,  this  is  consistent  with  the  eigenvalue 
representation  of  array  gain  where  small  eigenvalues  correspond  to  components  with  less 
noise. [17]  Since  Sx  must  be  positive  definite,  the  smallest  eigenvalue,  Amtn,  must  satisfy 

Amm  >  0  (2.41) 

Intuitively,  a  penalty  function  with  finite  (and  positive)  energy  throughout  the  scanning 
space  as  measured  at  the  sensor  is  specified  since  negative  or  zero  energy  systems  do  not 
possess  physical  interpretations  in  the  ocean  environment. 

With  this  interpretation  for  Sx,  positive  definite  matrices  can  easily  be  generated.  The 
value  of  Amtn  is  controlled  by  the  stability  factor.  By  generating  a  fixed  coherent  field  by 
equation  2.33,  Amtn  is  varied  by  adjusting  cr^  in  equation  2.32. 

While  the  requirement  of  equation  2.41  is  essential,  Amtn  must  also  be  considered 
jointly  with  Amax,  the  maximum  eigenvalue.  In  particular,  the  ratio  of  Xmax/Xmin  is  very 
important.  This  ratio  characterizes  the  dynamic  range  of  Sx.  A  large  dynamic  range 
places  more  demands  on  the  beamformer  than  a  small  dynamic  range;  therefore,  dynamic 
range  is  an  important  issue  when  considering  beamformer  performance.   Dynamic  range 
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is  applied  in  the  area  of  array  sensitivity  and  superdirectivity  and  is  considered  in  the 
next  section. 

2.4.2      Sensitivity  Ratio 

The  robustness  constraint  used  in  the  investigation  incorporates  the  principles  associ- 
ated with  the  white  noise  array  gain.  Sw,  related  to  Gw  by  equation  2.25,  is  a  measure  of 
sensitivity  to  uncorrelated  tolerance  errors.  [17]  To  ensure  robust  beamformer  performance 
a  constraint  is  applied  to  5W.[18]  The  robustness  constraint  is  applied  to  prevent  superdi- 
rective  performance.  Superdirective  performance  is  undesirable  because  a  superdirective 
array  is  very  sensitive  to  sensor  noise  and  uncorrelated  mismatch  errors.  Superdirectivity 
is  characterized  by  extremely  low  values  of  Gw. 

The  actual  constraint  applied  couples  Gw  to  a  parameter  called  the  array  sensitivity 
ratio  r)  which  is  defined  as 

n  =  JVWw  (2-42) 

Equation  2.42  is  derived  by  substituting  for  "Wmvdp  (as  defined  in  equation  2.19)  into 
equation  2.24.   After  some  matrix  manipulation  (since  Sx  is  Hermitian,  then  S"1  is  also 

Hermitian),  the  result  is 

_  letS-Mp 
GWimvdp  —  — r— o —  (.2.43) 

eT»x   e 

t  etSx2e  ,        , 

w'w  =  -, (2.44) 

(etS^e)2  l        } 


Calculating  wW  yields 


For  the  perfectly  matched  beamformer, 


w'w 


(2.45) 

(*w,MVDP 

Comparing  to  equation  2.25,  it  follows  that  the  magnitude  squared  of  the  element  weight- 
ing vector  is  a  measure  of  beamformer  sensitivity  (or  robustness.)  Since  Gw  is  always  less 
than  or  equal  to  N,  equation  2.42  is  the  ratio  of  the  uniformly  weighted  Gw  to  GWimvdp 
or  simply  a  measure  of  how  far  the  MVDP  beamformer  white  noise  gain  deviates  from 
N.  If  w  represents  a  uniformly  weighted  aperture,  then  rj  equals  1  and  the  beamformer  is 
robust.  As  the  beamformer  becomes  less  robust,  r/  increases  and  the  array  becomes  more 
sensitive. 
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The  sensitivity  ratio  as  defined  in  equation  2.42  is  a  very  convenient  measure  of  beam- 
former  sensitivity  since  w  is  available.  But  rj  must  be  related  to  superdirective  array 
performance.  Experience  shows  that  superdirectivity  occurs  at  sensitivity  ratios  in  the 
range  of  2.5  to  3.0;  this  leads  to  the  constraint 

r)  <  2.5  (2.46) 

The  value  of  77  is  controlled  by  varying  the  dynamic  range  of  Sr.  The  addition  of  a2N  to  the 
diagonal  elements  of  Sx  adds  a2N  to  each  eigenvalue  and  does  not  alter  the  eigenvectors. [18] 
Therefore,  the  inherent  structure  of  the  spectral  covariance  matrix  is  not  changed.  If 
^min  <<  Ij  then  adjusting  Amin  does  not  significantly  affect  Amaz  thus  controlling  the 
dynamic  range.  Large  dynamic  ranges  lead  to  superdirective  array  performance;  therefore, 
rj  is  controlled  by  limiting  the  dynamic  range  of  Sx. 

2.4.3     Directivity  Index 

Until  this  point,  the  analysis  of  the  array  performance  has  been  limited  to  an  ideal 
ocean;  one  that  has  a  perfectly  coherent  signal  in  incoherent  noise.  To  examine  the  effects 
in  a  more  realistic  ocean  environment,  the  requirement  of  incoherent  noise  is  relaxed  and 
the  case  of  isotropic  nose  is  considered. 

Directivity  index  is  defined  as  the  array  gain  of  a  perfectly  coherent  signal  in  isotropic 
noise.  Array  gain  can  also  be  defined  as  the  ratio,  in  decibel  units,  of  the  signal  to  noise 
ratio  of  an  array  to  the  signal  to  noise  ratio  of  a  single  element.  The  directivity  index  is 
obtained  by  assuming  a  perfectly  coherent  signal  which  has  cross  correlation  coefficients 
between  pairs  of  array  elements  of 

Mij  =  l       Vt,j  (2.47) 

and  isotropic  noise  with  cross  correlation  coefficients 

.     ,  sm(?f  \zi  -  zA) 

where  z,  is,  as  before,  the  position  of  the  ith  element.  Isotropic  noise  is  defined  as  having 
the  noise  power  per  unit  solid  angle  constant  in  all  directions.  Reference  [l]  describes  it  as 
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the  superposition  of  plane  waves  propagating  from  all  directions  with  uniform  statistical 
level.  Isotropic  noise  is  commonly  proposed  as  a  first  order  model  for  ambient  sea  noise. 
Under  these  assumptions,  the  directivity  index  (DI)  can  be  calculated  by 

^lOIog^fi^f'  (2.49) 

where  tu,-  is  the  weight  of  the  i**  element. [23]  DI  for  a  uniform  rectangular  array  with 
dimensions  corresponding  to  the  test  array  is  approximately  18  dB.[20] 

Use  of  DI  is  an  attempt  to  obtain  a  performance  measure  based  on  a  "more"  realistic 
ocean  than  previously  considered.  Since  "real"  ocean  noise  is  directional  in  both  the 
horizontal  and  vertical  and  ocean  multipaths  cause  reduced  signal  coherence,  the  use  of 
DI  has  practical  limitations.  However,  DI  is  commonly  used  and  is  useful  as  a  more 
realistic  measure  of  beamformer  performance  in  the  real  ocean. 
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Chapter  3 


Beamforming  Results 


3.1      Preliminaries 

Three  preliminary  issues  must  be  resolved  before  the  analysis  of  beam  pattern  designs 
can  proceed.  To  aid  in  beam  pattern  analysis,  a  visual  output  is  necessary.  The  relevant 
outputs  used  in  this  study  are  contour  plots  of  the  spatially  scanned  output,  the  beam 
pattern  and  the  square  magnitude  and  phase  of  the  element  weights.  The  contour  plots 
are  generated  from  a  uniform  grid  of  points.  For  the  array  output  and  beam  pattern 
contours,  the  grid  is  27  by  59.  For  the  element  weight  contours,  the  grid  is  10  by  20. 
Since  the  element  weights  are  complex,  it  is  convenient  to  plot  the  square  magnitude  of 
the  weights  and  a  corrected  phase.  The  square  magnitude  of  the  weights  are  scaled  to  a 
maximum  value  of  1  (or  0  dB).  The  corrected  phase  ($,-)  of  the  i*"  element  has  the  effect 
of  the  propagation  delay  across  the  array  removed.  This  is  calculated  by 

where  ^,-  is  the  element's  actual  phase,  ar  is  the  scan  direction,  and  zt-  is  the  element's 
position  vector.  All  four  plots  have  axes  corresponding  to  the  angular  position  in  the 
spatially  scanned  half  space.  For  the  magnitude  and  phase  plots,  these  axes  relate  the 
element  positions  in  terms  of  their  location  in  the  scanning  field. 

The  next  issue  is  the  verification  of  the  beamforming  algorithm.  The  test  used  to 
verify  proper  operation  is  a  plane  wave  incident  at  0°  elevation  and  90°  azimuth  (true 
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broadside  case)  in  the  presence  of  sensor  noise  only.  For  this  test  the  expected  results 
are  a  sensitivity  ratio  of  1  (i.  e.  GWimvdp  =  N  =  200)  and  a  uniform  weighting  vector. 
The  MVDP  beamformer  results  are  as  expected  with  DI  =  39.2.  The  array  output  for 
this  plane  wave  case  is  shown  in  figure  3-1.  Observe  that  the  array  output  correctly 
detects  the  plane  wave.  Figure  3-1  is  also  useful  in  that  it  provides  some  insight  into  the 
azimuthal  and  elevational  resolution  of  the  array  processor.  The  beam  pattern  resulting 
from  the  plane  wave  input  is  the  conventional  beam  pattern  steered  to  0°  elevation  and  90° 
azimuth;  this  is  shown  in  figure  3-2.  In  figure  3-2,  observe  that  there  is  unity  gain  in  the 
scan  direction  as  required  in  the  MVDP  derivation.  Also,  observe  the  first  sidelobes  (in 
both  elevation  and  azimuth)  are  at  -13  dB  which  is  as  expected  for  a  uniformly  weighted 
rectangular  array.  The  second  sidelobes  in  azimuth  at  -18  dB  also  correlate  well  with  the 
rectangular  array  approximation.  These  two  plots  are  useful  as  points  of  comparison  for 
the  broadside  results.  Results  are  also  included  for  steered  patterns,  in  particular,  patterns 
steered  to  0°  elevation  and  102°  azimuth.  For  comparison  purposes  in  anticipation  of  the 
steered  results,  the  conventional  beam  pattern  steered  to  0°  elevation  and  102°  azimuth 
is  included  in  figure  3-3. 

The  last  issue  concerns  the  sufficiency  of  1/4  BW  scanning.  When  a  simple  penalty 
function,  set  up  to  ensure  directions  between  scanning  increments  are  present,  is  input  to 
the  beamformer  for  the  cases  of  1/8  BW ,  1/4  BW  and  1/2  BW  scanning,  no  observable 
effects  due  to  incident  plane  waves  between  scan  directions  are  present.  Visually,  the 
finer  the  scanning  increment,  the  better  (qualitatively)  the  output  appears  which  is  as 
expected.  These  results  confirm  the  sufficiency  of  ^BW  scanning. 

3.2      Broadside  Results 

The  analysis  of  this  investigation  centers  on  a  systematic  evaluation  of  the  effects  of 
different  penalty  functions  on  beamformer  performance  (as  measured  by  the  sensitivity 
ratio  and  directivity  index)  and  beam  pattern  structure.  The  parameters  specifying  the 
penalty  function  are  listed  in  table  3.1  and  are  identified  in  chapter  2.3.  The  center  region 
widths  are  easily  defined;  they  are  simply  the  widths  of  the  plateau.  The  transition  region 
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Figure  3-1:  array  output  for  plane  wave  case 
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Figure  3-2:  conventional  beam  pattern  steered  to  0°  elevation  and  90°  azimuth 
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Figure  3-3:  conventional  beam  pattern  steered  to  0°  elevation  and  102°  azimuth 
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parameter 

symbol 

minimum  eigenvalue 

^mm 

center  region  width 

elevation 

A0C 

azimuth 

A<f>c 

transition  region  width 

elevation 

A$t 

azimuth 

Mt 

Table  3.1:  Penalty  Function  Parameters 

widths  are  more  difficult  to  define.  The  widths  are  controlled  by  varying  a\  and  a\  in 
equation  2.38.  From  this,  A9t  and  A<f>t  are  defined  in  terms  of  the  3  dB  down  points 
along  the  major  axes.  All  region  widths  are  measured  in  degrees.  The  investigation 
concentrates  on  penalty  functions  centered  about  0°  elevation  and  90°  azimuth;  the  bulk 
of  the  results  come  from  the  study  of  the  broadside  case.  This  scan  direction  is  used  for 
several  reasons.  First,  patterns  produced  are  symmetric  and  are  very  easy  to  analyze. 
A  shift  off  of  this  scan  direction  produces  very  asymmetric  patterns;  these  patterns  are 
difficult  to  analyze.  Second,  this  scan  direction  produces  patterns  that  can  be  intuitively 
compared  to  a  uniform  rectangular  array  steered  to  broadside.  Finally,  it  is  easier  to 
control  the  array  output  for  method  B  penalty  functions  with  the  symmetric  patterns. 
Since  the  two  methods  investigated  have  different  behavior  and  characteristics,  they  will 
be  addressed  separately. 

3.2.1      Method  A 

The  first  parameter  to  be  resolved  for  both  methods  is  the  selection  of  Amin  or  more 
specifically,  the  selection  of  the  dynamic  range  of  PF.  From  section  2.4.1,  the  larger 
Amtn  can  be  correlated  with  larger  values  of  a2N.  The  more  the  a2N  term  dominates  in 
equation  2.32,  the  closer  the  penalty  function  resembles  a  sensor  noise  only  case.  This 
leads  to  the  conventional  beam  pattern.  Intuitively,  the  most  interesting  cases  are  for 
small  Xmin  or  a  large  dynamic  range.  Experience  gained  from  working  with  method  A 
penalty  functions  shows  that  this  method  has  very  predictable  performance  with  respect 
to  dynamic  range.  For  the  penalty  function  parameters  investigated  with  method  A,  the 
dynamic  range  which  leads  to  a  sensitivity  ratio  of  approximately  2  is  on  the  order  of  104 
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minimum  eigenvalue  (Amin) 

dynamic  range 

1-1CT3 

1.1-  104 

2-1CT3 

5.6  -103 

3  •  1(T3 

3.7  •  103 

5  •  icr3 

2.2  •  103 

Table  3.2:  Minimum  Eigenvalues  and  Associated  Dynamic  Ranges  Tested  for  Method  A 
Penalty  Functions 

with  a  corresponding  Amin  on  the  order  of  10-3.   For  purposes  of  comparison,  Amtn  for 
method  A  is  restricted  to  the  values  listed  in  table  3.2. 

With  a  dynamic  range  established,  15  cases  are  used  to  explore  the  effects  of  transition 
and  center  region  widths  on  method  A  penalty  functions.  These  15  cases  are  summarized 
in  table  3.3.  Cases  1A  through  5A  examine  the  effects  of  center  region  width  in  azimuth. 
The  best  results  in  terms  of  DI  and  first  sidelobe  level  occur  with  case  2A.  For  this  case  the 
first  sidelobe  level  in  azimuth  is  approximately  -24.9  dB  and  in  elevation  is  approximately 
-22.8  dB  with  DI  =  37.9  and  tj  =  2.09.  For  all  cases  1A  through  5A,  changing  A<f>c 
has  no  effect  on  the  first  sidelobe  level  in  elevation.  Cases  6A  through  HA  examine  the 
effects  of  center  region  width  in  elevation.  Widths  beyond  52.5°  cause  the  first  sidelobe 
to  disappear.  This  results  in  wider  main  beam  widths;  therefore,  cases  beyond  52.5°  are 
not  considered.  The  best  results  occur  with  case  HA.  The  resulting  first  sidelobe  level 
in  azimuth  is  approximately  -32.9  dB  and  in  elevation  is  approximately  -36.2  dB  with 
DI  =  37.4  and  r\  =  1.62.  The  eigenvalue  performance  characteristics  for  case  HA  in 
terms  of  77  and  DI  are  shown  in  figure  3-4.  As  opposed  to  cases  1A  through  5A,  changing 
A9C  significantly  affects  first  sidelobe  levels  in  azimuth  where  up  to  8  dB  differences  are 
observed.  This  implies  a  partial  coupling  between  azimuth  and  elevation  first  sidelobe 
levels.  Using  case  11 A  as  the  comparison  point,  cases  12 A  and  13 A  examine  the  effects 
of  transition  region  width  in  azimuth  while  cases  14A  and  15A  examine  the  effects  of 
transition  region  width  in  elevation.  Transition  region  width  appears  to  serve  as  a  fine 
tuner  in  defining  an  effective  method  A  penalty  function  width.  The  best  results  occur 
for  case  11  A;  therefore,  the  array  output,  beam  pattern,  square  magnitude  and  corrected 
phase  for  the  elements  weights  are  shown  in  figures  3-5  through  3-8.  Figure  3-5  is  the  array 
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case 

A9l 

Ml 

A0j 

A<rf 

1A 

15 

15 

16.2 

16.2 

2A 

15 

22.5 

16.2 

16.2 

3A 

15 

30 

16.2 

16.2 

4A 

15 

37.5 

16.2 

16.2 

5A 

15 

45 

16.2 

16.2 

6A 

7.5 

22.5 

16.2 

16.2 

7A 

22.5 

22.5 

16.2 

16.2 

8A 

30 

22.5 

16.2 

16.2 

9A 

37.5 

22.5 

16.2 

16.2 

10A 

45 

22.5 

16.2 

16.2 

11A 

52.5 

22.5 

16.2 

16.2 

12A 

52.5 

22.5 

16.2 

11.7 

13A 

52.5 

22.5 

16.2 

21.0 

14A 

52.5 

22.5 

8.8 

16.2 

15A 

52.5 

22.5 

28.6 

16.2 

*all  widths  in  degrees 


Table  3.3:  Method  A  Case  Study  Summary 

processor's  representation  of  the  penalty  function.  Figure  3-6  has  significantly  reduced 
sidelobes  compared  to  the  conventional  beam  pattern  shown  in  figure  3-2.  The  weights, 
as  defined  in  figures  3-7  and  3-8,  do  not  possess  a  structure  that  is  readily  identified  with 
a  conventional  weighting. 

Based  on  the  method  A  performance  results,  the  penalty  function  design  process  for 
method  A  is  summarized  by  the  following  rules: 
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Figure  3-4:  Eigenvalue  Performance  Results  for  Case  11 A 
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Figure  3-5:  array  output  for  case  11A,  Amtn  =  10 
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Figure  3-6:  beam  pattern  for  case  11  A,  Amtn  =  10 
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Figure  3-7:  square  magnitude  of  element  weights,  case  11A,  Amtri  =  10 
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Figure  3-8:  corrected  phase  of  element  weights,  case  11A,  Amtn  =  10  3 
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•  Selection  of  center  region  width  is  the  most  critical  parameter  choice  and  should  be 
on  the  order  of  the  width  of  the  main  beam  of  the  conventional  beam  pattern. 

•  Center  region  widths  in  directions  with  a  large  aperture  should  be  selected  first  since 
these  widths  appear  to  be  uncoupled  with  directions  with  short  apertures. 

•  The  transition  region  widths  should  be  selected  to  fine  tune  the  effective  penalty 
function  width  to  yield  a  smooth  transition  over  a  region  less  than  the  width  of  the 
center  region. 

•  Select  a  dynamic  range  that  prevents  superdirective  performance.  Large  dynamic 
ranges  (but  still  not  superdirective)  usually  result  in  lower  first  sidelobe  levels,  wider 
main  beams,  and  smaller  DI,  so  some  trade  off  may  be  desired  in  the  selection  of 
dynamic  range. 

Using  these  guidelines,  a  significant  reduction  in  sidelobe  levels  can  be  achieved  using  the 
method  A  penalty  function  design  process. 

3.2.2      Method  B 

As  was  the  case  for  method  A,  the  dynamic  range  for  method  B  penalty  functions  must 
first  be  established.  This  issue  is  not  as  easily  addressed  for  method  B  as  it  is  for  method 
A.  Although  the  most  interesting  cases  are  still  for  small  Amtfl,  method  B  penalty  functions 
exhibit  unpredictable  behavior  in  such  a  way  that  dynamic  range  must  be  considered  as 
a  separate  issue  for  each  case.  Therefore  to  cover  the  majority  of  possibilities,  the  values 
of  Amm  listed  in  table  3.4  are  considered.  These  values  of  Amtn  correspond  to  dynamic 
ranges  varying  from  106  to  101. 

The  effects  of  transition  and  center  region  widths  for  method  B  penalty  functions 
are  explored  in  13  cases  which  are  summarized  in  table  3.5.  Cases  4B  and  5B  have 
performance  characteristics  in  terms  of  rj  and  DI  versus  Amtn  similar  to  method  A;  all 
are  superdirective  at  large  values  of  dynamic  range.  The  remaining  cases  have  different 
performance  characteristics.  These  cases  can  not  be  driven  superdirective  even  for  very 
large  dynamic  ranges.  In  fact,  these  cases  exhibit  a  stable  performance  region  where  no 
gain  in  sidelobe  reduction  is  achieved  with  increasing  dynamic  range.  The  performance 
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minimum  eigenvalue  (Amtn) 

dynamic  range 

l(Tb 

9.5  •  105 

1(T5 

9.5  •  104 

10-4 

9.5  •  103 

1CT3 

9.5  •  10'2 

KT2 

9.5  •  101 

ict1 

9.5  •  10° 

Table  3.4:  Minimum  Eigenvalues  and  Associated  Dynamic  Ranges  Tested  for  Method  B 
Penalty  Functions 

results  for  case  9B  (this  case  displays  the  best  results  of  all  13  cases)  are  shown  in  figure 
3-9.  Observe  the  stable  operating  region  for  Amtn  <  10-3. 

Cases  IB  through  5B  examine  the  effects  of  A<pc.  For  cases  4B  and  5B,  A<f>c  is  too 
wide  causing  the  first  sidelobes  to  be  emphasized  in  addition  to  the  main  lobe.  This  is 
the  cause  of  the  superdirective  performance  results.  Case  3B  results  in  the  best  sidelobe 
control  with  first  sidelobe  level  in  azimuth  equal  to  -49.7  dB,  first  sidelobe  level  in  elevation 
equal  to  -17.3  dB,  Amtn  =  10-4,  77  =  1.6,  and  DI  =  34.6  dB.  No  effects  are  observed  on 
sidelobe  levels  in  elevation  for  cases  IB  through  5B. 

Cases  6B  through  8B  examine  the  effects  of  A0C.  All  cases  produce  first  sidelobe  levels 
within  1  dB  of  each  other,  the  best  being  case  7B  with  a  first  sidelobe  level  in  elevation  of 
-18  dB.  Sidelobe  levels  in  elevation  are  generally  insensitive  to  method  B  penalty  function 
parameter  adjustment. 
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Figure  3-9:  Eigenvalue  Performance  Results  for  Case  9B 
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case 

A0* 

Ml 

A0* 

A<p\ 

IB 

15 

15 

16.2 

16.2 

2B 

15 

22.5 

16.2 

16.2 

3B 

15 

30 

16.2 

16.2 

4B 

15 

37.5 

16.2 

16.2 

5B 

15 

45 

16.2 

16.2 

6B 

7.5 

22.5 

16.2 

16.2 

7B 

22.5 

22.5 

16.2 

16.2 

8B 

30 

22.5 

16.2 

16.2 

9B 

22.5 

30 

16.2 

16.2 

10B 

22.5 

30 

16.2 

11.7 

11B 

22.5 

30 

16.2 

21.0 

12B 

22.5 

30 

8.8 

16.2 

13B 

22.5 

30 

28.6 

16.2 

*all  widths  in  degrees 


Table  3.5:  Method  B  Case  Study  Summary 

Case  9B  represents  the  best  results  of  the  center  region  width  study  and  is  used  as  the 
comparison  case  for  the  transition  region  width  study.  Cases  10B  through  13B  investigate 
the  transition  region  width  issue.  Transition  region  effects  for  method  B  are  similar  to 
those  effects  for  method  A;  transition  region  width  serves  as  fine  tuning  for  the  effective 
penalty  function  width. 

Case  9B  represents  the  best  results  for  all  15  cases.  The  array  output,  beam  pattern, 
square  magnitude  and  corrected  phase  of  the  element  weights  are  shown  in  figures  3-10 
through  3-13.  Figure  3-10  shows  the  array  processor's  interpretation  of  the  method  B 
penalty  function.  The  beam  pattern  in  figure  3-11  has  significantly  reduced  sidelobes  in 
azimuth  compared  to  the  conventional  beam  pattern  found  in  figure  3-2,  but  has  only 
minor  improvements  in  sidelobe  levels  in  elevation.  The  element  weights  shown  in  figures 
3-12  and  3-13  have  a  structure  similar  to  conventional  weights:  high  in  the  middle  and 
tapering  towards  the  ends. 

The  method  B  investigation  also  revealed  several  important  issues  characteristic  of  the 
method  B  penalty  function.  One,  the  beamformer  tends  to  emphasize  edges;  in  particular, 
the  edges  corresponding  to  the  transition  region.  If  the  center  region  is  wide  enough,  the 
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Figure  3-10:  array  output,  case  9B,  Amtn  =  10~4 
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Figure  3-12:  square  magnitude  of  element  weights,  case  9B,  Amin  =  10~4 
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Figure  3-13:  corrected  phase  of  element  weights,  case  9B,  Am,n  =  10~4 
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edges  are  located  in  the  vicinity  of  the  first  sidelobe  of  the  conventional  beam  pattern. 
This  is  shown  in  figure  3-14,  the  case  5B  array  output.  This  results  in  the  superdirective 
beam  pattern  found  in  figure  3-15.  More  significantly,  the  element  weights  appear  to  have 
a  threshold  where  the  structure  changes  suddenly  to  a  significantly  different  shape.  This 
is  shown  in  figures  3-16  and  3-17. 

The  penalty  function  design  process  for  method  B  is  summarized  by  the  following 
rules. 

•  Select  a  center  region  width  on  the  order  of  the  width  of  the  main  beam  of  the 
conventional  beam  pattern. 

•  Selection  of  the  center  region  width  involves  a  trade  off  against  the  array  aperture. 
For  those  directions  with  large  apertures  (azimuth),  the  full  main  beam  width  should 
be  used.  For  those  directions  with  short  apertures  (elevation),  the  center  region 
width  must  be  limited  to  something  less  than  the  main  beam  width  due  to  the  wide 
main  beam  and  low  sensitivity  to  parameter  adjustment. 

•  Select  transition  region  widths  to  fine  tune  the  effective  penalty  function  width  and 
to  yield  a  smooth  transition  over  a  region  less  than  the  width  of  the  center  region. 

•  Select  the  smallest  possible  dynamic  range  that  is  located  in  a  stable  operating  region 
if  a  stable  region  exists.  Otherwise,  choose  the  largest  dynamic  range  consistent  with 
superdirectivity  and  DI  considerations. 

Using  these  guidelines  for  method  B  designs,  a  significant  reduction  in  levels  can  be 
achieved  with  the  more  notable  results  occurring  in  directions  of  large  array  aperture. 

3.3      Steered  Results 

The  problem  encountered  with  steered  penalty  functions  is  the  introduction  of  asym- 
metry into  the  problem.  This  asymmetry  has  very  little  effect  on  method  A  and  has 
very  significant  effects  on  method  B.  Only  limited  results  are  presented  here;  a  single  case 
steered  to  0°  elevation  and  102°  azimuth  is  illustrated.  The  intent  is  to  only  provide  an 
intuitive  feel  for  steered  penalty  function  performance. 

59 


scale  in  dB 


■>*     -0;     t^     O    CO     CO    0>    C\J 

iriirii^ocji-cooi    zi 
NWNnonnn    f 


"*   r- 


q  n  id   01  w  * 

B    Q     N    U)    N 

■<t   Co    in   in   in 


a    r*«   *j    k>i 

25    Cvi    CO     CO 


CO     ID    Ol    CNI 

rt  rt  n   * 


*son«joi«« 


ui    in   in 


90 


70 


50 


30 


CO 

(1) 

CD 

10 

D) 

© 

■D 

C 

-10 

o 

*-• 

CO 

> 

CD 
CD 

-30 

■50 


■70 


-90 


60         80        100       120       140 
azimuth  (degrees) 


160      180 


Figure  3-14:  array  output,  case  5B,  Amtn  =  10 
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Figure  3-16:  square  magnitude  of  element  weights,  case  5B,  Amtn  =  10" 


62 


scale  in  pi  radians 

."£    ."£    S2     <M    £     t~~    °    CO     ID    Ol    N    ■*    s 

12  ~  ^  01  ^  n  ;;  ra  Ln<Nor~Tj- 
pop  o  ■:  w  ^  o  t  ui  10  (o  s 
o   o   o    ©   o   <->   o   <  ' 


9  9 


o   n   co 

CM    O     ID 
CO     CO     Ol 


9   9  9 


cd  en  cm  "J  r- 
>  t-  en  <c  co 
°    p   p    - 


^j    «   w 

■"".     CM     W     CO 

o  o  o  p  o 


CO     CO      CO     CM 

CO     Ul      IM      " 


If)     ID     CO     N 


O     CO     CO 
CM     O)     CO 

co    co   o> 


9  9  9  9   9  9 


mmnn 


0    20    40    60    80   100   120   140   160 

azimuth  (degrees) 


180 


Figure  3-17:  corrected  phase  of  element  weights,  case  5B,  Amin  =  10-4 
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3.3.1  Method  A 

Virtually  no  effects  due  to  steers  are  encountered.  When  the  case  9A  penalty  function 
is  steered  to  0°  elevation  and  102°  azimuth,  no  distortion  of  the  penalty  function  in  terms 
of  the  array  output  is  observed.  This  is  shown  in  figure  3-18.  The  really  nice  feature  of  the 
method  A  procedure  is  that  steering  the  penalty  function  has  essentially  no  effect  on  the 
sidelobe  levels  compared  to  the  broadside  case.  This  is  shown  in  figure  3-19.  Again,  this 
beam  pattern  has  significantly  better  sidelobe  levels  than  the  steered  conventional  beam 
pattern  found  in  figure  3-3.  The  square  magnitude  and  corrected  phase  of  the  element 
weights  for  the  steered  case  9A  are  shown  in  figures  3-20  and  3-21.  Again,  the  aperture 
weights  are  not  characteristic  of  any  traditional  shading. 

3.3.2  Method  B 

The  effects  of  asymmetry  on  the  method  B  penalty  function  are  manifested  in  the 
concept  of  edge  detections.  The  asymmetry  causes  the  edges  in  directions  of  shorter 
apertures  to  be  more  emphasized  which  causes  a  significant  sidelobe  to  be  placed  in  the 
location  of  the  short  aperture  edge.  This  is  observed  when  case  9B  is  steered  to  0°  elevation 
and  102°  azimuth;  the  array  is  superdirective  in  this  case.  As  the  penalty  function  center 
region  is  made  wider  and  as  the  array  is  steered  further  off  broadside,  the  asymmetry 
effects  become  more  pronounced.  (If  the  center  region  is  too  wide,  a  threshold  effect 
is  again  observed  where  the  weights  suddenly  change  character  into  a  dumbbell  shape.) 
The  placement  of  significant  sidelobes  at  the  locations  of  edges  is  similar  to  that  which  is 
observed  for  the  broadside  case  5B.  This  effect  is  illustrated  by  observing  the  array  output 
and  beam  patterns  for  case  5B  in  figures  3-14  and  3-15  where  both  edges  are  emphasized 
instead  of  just  one  edge  in  the  direction  of  short  aperture  as  is  found  in  the  steered  case 
9B. 

Based  on  these  observations,  the  method  B  penalty  function  procedure  for  steered 
patterns  must  be  modified.  The  key  is  to  reduce  the  center  region  width  to  a  point  where 
the  beamformer  cannot  distinguish  one  edge  from  the  other;  rather  than  two  nulls  in 
the  array  output,  there  appears  just  one  center  null.  If  A<f>c  is  reduced  to  0°,  then  more 
reasonable  results  are  obtained.  The  array  output,  beam  pattern  and  square  mag  n  nude 
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Figure  3-18:  array  output,  steered  method  A  penalty  function 
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Figure  3-19:  beam  pattern,  steered  method  A  penalty  function 
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Figure  3-20:  square  magnitude  of  element  weights,  steered  method  A  penalty  function 
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Figure  3-21:  corrected  phase  of  element  weights,  steered  method  A  penalty  function 


68 


and  corrected  phase  of  the  element  weights  for  the  zero  center  region  width  steered  to 
0°  elevation  and  102°  azimuth  are  shown  in  figures  3-22  through  3-25.  Although  these 
results  are  not  as  good  as  the  broadside  case,  the  beam  pattern  still  has  reduced  sidelobes 
compared  to  the  steered  conventional  beam  pattern.  Better  performance  may  be  obtained 
if  some  small  center  region  width  is  inserted. 
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Figure  3-22:  array  output,  steered  method  B  penalty  function 
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Figure  3-23:  beam  pattern,  steered  method  B  penalty  function 
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Figure  3-24:  square  magnitude  of  element  weights,  steered  method  B  penalty  function 
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Figure  3-25:  corrected  phase  of  element  weights,  steered  method  B  penalty  function 
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Chapter  4 


Summary  and  Conclusions 


This  investigation  features  a  systematic  procedure  to  design  beam  patterns  for  a 
general  multidimensional  irregular  array,  an  area  where  little  previous  work  has  been 
accomplished.  The  procedure  uses  a  "penalty  function"  input  to  an  MVDP  beamformer. 
The  proper  penalty  function  penalizes  high  sidelobes  and  encourages  the  main  lobe.  Two 
different  methods  of  penalty  function  design  are  investigated.  Method  A  achieves  sidelobe 
reduction  by  placing  nulls  throughout  the  scanning  space  encompassing  the  entire  sidelobe 
region.  This  method  is  characterized  by  fully  defining  the  scanning  region;  this  leads 
to  very  predictable  and  stable  performance.  Method  B  achieves  sidelobe  reduction  by 
emphasizing  the  main  lobe.  This  is  achieved  by  placing  anti-nulls  in  the  main  beam 
region.  Method  B  is  characterized  by  a  sparsely  defined  scanning  region;  this  leads  to 
sometimes  unpredictable  and  unstable  performance.  The  design  procedure  is  a  four  step 
process: 

1.  Determine  the  main  lobe  beam  width  for  the  conventional  beam  pattern.  The  main 
lobe  width  is  best  defined  in  terms  of  the  angular  separation  between  the  first  nulls. 

2.  Select  the  center  region  width  based  on  the  main  lobe  beam  width  and  the  design 
method  (A  or  B). 

3.  Fine  tune  the  effective  penalty  function  width  by  selecting  an  appropriate  transition 
region  width. 
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Figure  4-1:  Plot  of  First  Sidelobe  versus  Center  Region  Width,  Method  A.  Conventional 
main  lobe  widths  in  terms  of  angular  separation  between  first  nulls  are  approximately  24° 
in  azimuth  and  50°  in  elevation. 

4.  Stabilize  the  beamformer  by  the  addition  of  a  sensor  noise  term  in  such  a  way  to 
prevent  superdirective  performance. 

Selection  of  the  center  region  width  is  the  critical  step.  There  is  an  optimum  center  region 
width  with  respect  to  the  resultant  first  sidelobe  level.  This  is  illustrated  for  both  methods 
as  applied  to  the  array  used  throughout  the  investigation;  these  results  are  shown  in  figures 
4-1  and  4-2.  The  figures  show  that  a  minimum  in  sidelobe  level  occurs  for  each  direction 
except  in  elevation  for  method  A  where  increasing  the  center  region  width  eliminates  the 
first  sidelobe. 

Methods  A  and  B  each  have  distinctive  characteristics.  The  character  of  each  method 
must  be  evaluated  when  considering  which  method  to  apply. 

•  Method  A  has  roughly  equivalent  sidelobe  performance  in  elevation  and  azimuth. 
Method  B  performance  in  azimuth  is  significantly  better  than  performance  in  ele- 
vation. 

•  In  terms  of  first  sidelobe  level  overall  performance  in  elevation  is  better  with  method 
A;  overall  performance  in  azimuth  is  better  with  method  B. 

•  Method  A  results  in  a  narrower  main  lobe  compared  to  method  B  which  results  in 
better  DI  performance  for  method  A. 
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Figure  4-2:  Plot  of  First  Sidelobe  versus  Center  Region  Width,  Method  B.  Conventional 
main  lobe  widths  in  terms  of  angular  separation  between  first  nulls  are  approximately  24° 
in  azimuth  and  50°  in  elevation. 

•  Method  A  superdirective  performance  is  characterized  by  high  sidelobes  in  the  un- 
observed k  space.  Method  B  superdirective  performance  is  characterized  by  high 
sidelobes  in  the  observable  k  space. 

•  Method  A  is  easily  steered;  the  method  B  procedure  must  be  modified  when  steered. 

•  For  method  A,  azimuth  performance  is  partially  coupled  to  elevation  performance. 
For  method  B,  azimuth  performance  and  elevation  performance  are  uncoupled. 

•  Method  B  has  stable  operating  regions  where  performance  is  independent  of  penalty 
function  dynamic  range. 

•  The  method  B  procedure  has  a  threshold  effect  where  the  aperture  weights  suddenly 
change  character  from  shapes  with  one  global  maximum  to  shapes  with  two  or  more 
local  maximum. 

The  investigation  is  notably  sparse  of  results  for  steered  arrays;  however,  with  the 
use  of  some  simple  cases,  significant  sidelobe  reduction  is  achieved  for  the  examined  case. 
Only  the  effects  of  steers  in  azimuth  are  investigated  because  of  the  insensitivity  found  in 
elevation.  Considering  the  case  presented,  the  general  requirements  for  steered  penalty 
functions  are  illustrated. 
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Additional  investigation  into  several  aspects  not  examined  are  warranted.  Method  A 
is  the  design  method  of  choice,  however,  it  would  be  beneficial  to  obtain  sidelobe  reduction 
in  azimuth  comparable  to  that  found  in  method  B.  One  proposal  to  achieve  these  levels 
would  be  the  incorporation  of  a  two  level  background  region  which  more  heavily  penalizes 
sidelobes  in  azimuth. 

Finally,  one  major  issue  which  has  been  only  surficially  addressed  concerns  the  robust- 
ness or  sensitivity  of  the  design  procedure  in  terms  of  the  uncertainty  in  the  amplitudes 
and  phases  of  the  response.  This  issue  is  not  well  understood  and  yet  is  so  important 
when  characterizing  the  performance  of  multidimensional  arrays.  Further  work  on  this 
robustness  issue  is  required  and  is  the  next  step  in  refining  the  "penalty  function"  design 
procedure. 
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Appendix  A 

Spectral  Covariance  Matrix 
Inversion 


When  considering  an  MLM  beamforming  approach,  the  problem  of  solving  a  large 
system  of  equations  quickly  arises.  At  the  outset  of  this  undertaking,  it  was  hoped  ex- 
ploiting the  embedded  Toeplitz  structure  of  the  spectral  covariance  matrix  would  lead  to 
significant  savings  in  computational  time.  What  follows  are  the  results  of  the  investigation 
into  Toeplitz  matrix  inversion  techniques. 
Consider  the  matrix  problem 

Ax  =  b  (A.l) 

where  A  and  b  are  known.  The  standard  matrix  solution  method  which  solves  any  nonsin- 
gular  system  is  Gaussian  elimination. [24]  The  method  is  well  documented  and  numerous 
algorithms  are  readily  available  which  efficiently  incorporate  the  method,  particularly  the 
LINPACK  routines. [12]  The  problem  with  this  method  is  that  it  is  generally  the  slowest 
of  the  standard  methods.  For  an  n  by  n  real  matrix,  the  number  of  multiplications  for  the 
Gaussian  elimination  method  is  proportional  to  n3.  A  Toeplitz  inversion  method  using 
a  more  efficient  bordering  method  was  developed  by  Levinson  with  the  number  of  multi- 
plications proportional  to  n2.[8]  The  goal  was  to  incorporate  these  savings  into  the  more 
complicated  structure  of  the  spectral  covariance  matrix. 

The  multidimensional  beamforming  problem  produces  a  more  complicated  matrix 
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scanning  pattern 

method 

time(sec) 

aggregate  pattern 

Gaussian  elimination 

4400 

Toeplitz 

4950 

single  scan 

Gaussian  elimination 

265 

Toeplitz 

160 

Table  A.l:  Computational  Times  for  Matrix  Solution  Methods 

structure  as  discussed  in  chapter  2.2.2.  Two  level  Toeplitz  matrix  solution  methods  suit- 
able for  this  complicated  matrix  structure  have  been  developed. [25, 13]  The  method  used 
is  from  the  Toeplitz  Package  Users  Guide. [13]  The  number  of  multiplications  for  a  real 
matrix  with  this  method  is  approximately  2M3L2.  For  the  particular  array  structure, 
M  —  10  and  L  =  20.  This  represents  a  factor  of  10  savings  from  the  Gaussian  elimination 
method  with  n  =  200.  Even  with  the  overhead  resulting  from  the  permutations  required 
to  obtain  the  proper  structure  for  input  into  the  Toeplitz  algorithm,  a  significant  savings 
is  expected.  Two  factors  have  not  been  accounted  for  up  to  this  point.  One,  the  matrices 
are  complex.  Two,  the  system  of  equations  needs  to  be  solved  1593  times. 

The  two  methods  were  tested  against  each  other  with  a  1600  scanning  direction  pat- 
tern, where  1600  vice  1593  scanning  directions  were  used  to  accommodate  permutations 
of  the  right  hand  side  of  equation  A.l,  consistent  with  the  permutation  of  A.  The  tim- 
ing results  were  obtained  on  a  Digital  Micro vax  II  computer  system.  The  run  times  for 
all  scanning  directions  and  for  just  one  scanning  direction  are  summarized  in  table  A.l. 
Although  the  actual  beamforming  was  performed  on  an  Alliant  FX40  computer  system, 
the  time  performances  for  both  methods  were  similar.1 

The  superior  performance  of  the  Gaussian  elimination  method  for  the  full  problem  can 
be  attributed  to  the  LU  decompositon  method  used.  Considering  the  single  run  times, 
the  full  Toeplitz  procedure  is  dependent  on  b  and  must  be  repeated  for  each  new  right 
hand  side  vector.  Only  a  portion  of  the  Gaussian  elimination  procedure  is  dependent  on 
b.  The  matrix  decomposition  need  only  be  accomplished  once  while  the  back  substitution 
using  the  different  right  hand  side  vectors  can  be  applied  to  the  same  decomposition. 


'No  individual  run  times  were  examined  on  the  Alliant,  only  aggregate  times.  The  Alliant  times  were 
approximately  65  times  faster  than  the  Microvax. 
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In  conclusion,  the  Toeplitz  method  is  more  efficient  when  the  matrix  problem  only 
requires  a  few  right  hand  side  vectors;  Gaussian  elimination  is  more  efficient  when  a  large 
number  of  right  hand  sides  must  be  evaluated. 
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Appendix  B 


Dolph-Chebyshev  Weights 


Aperture  weights  have  been  presented  for  method  A  and  method  B  penalty  functions. 
As  previously  noted,  method  A  aperture  weights  have  no  apparent  special  structure,  but 
method  B  aperture  weights  have  a  structure  which  may  be  compared  to  a  conventional 
weight.  Since  the  MVDP  beamformer  is  optimum  in  the  minimum  variance  sense,  a 
useful  comparison  might  be  made  with  Dolph-Chebyshev  weights  which  are  optimum  in 
the  relationship  between  minor  lobe  level  and  main  lobe  width.  Using  the  procedures 
outlined  in  reference  [26],  the  aperture  weighting  for  a  uniform  rectangular  array  can  be 
calculated  with  input  parameters  of  side-lobe  height  in  both  azimuth  and  elevation.  The 
results  for  -60  dB  side-lobes  in  azimuth  and  -22  dB  side-lobes  in  elevation  are  shown  in 
figure  B-l.  The  Dolph-Chebyshev  aperture  weights  in  figure  B-l  are  remarkably  similar 
to  the  aperture  weights  for  case  9B  in  figure  3-12,  particularly  away  from  the  edges  of 
the  array.  Furthermore,  the  side  lobe  levels  for  the  case  9B  beam  pattern  in  figure  3- 
11  can  be  characterized  by  -60  dB  in  azimuth  and  -22  dB  in  elevation.  Based  on  these 
results,  method  B  penalty  function  aperture  weights  can  be  favorably  compared  to  Dolph- 
Chebyshev  aperture  weights. 
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Figure  B-l:  Dolph-Chebyshev  square  magnitude  element  weights 
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